Because of the difference of impedance of each medium involved in the MIM structure, the waves that propagate inside this structure will experience multiple reflections and transmissions at each interface. These reflections can occur an infinite number of times, and their summation is represented by a geometric series. Transmission “loops” can be defined inside the structure in order to enumerate every possible path of light inside the structure. In a multilayer system like the MIM structure, each layer and combination of layers will behave as a possible “loop” of the transmission. represents the different loops possible in the four-interface MIM structure.

Figure 3: (A) The different loops inside the MIM structure. (B) The expansion of loop *f*_{23} by taking into account the first diffraction orders.

The total transmission becomes a complex summation of all possible paths, which will define how the fields inside each layer interact with each other. Light can be transmitted directly, as represented by the *t*^{d} arrow in , or can go through a series of loops of the system before being transmitted. Loops can also exist inside other loops.

Waves inside the slits can propagate as different modes of the parallel plate waveguide, so all possible loops containing these regions can be expanded into the different combinations of the guided modes considered. Furthermore, plan waves in region 3 can propagate as diffraction orders of the periodic system, as described by the indices of the reflection and transmission coefficients ${\tau}_{\alpha \mathrm{,}k}^{34},$${\tau}_{\alpha \mathrm{,}k}^{23},$${\rho}_{{k}_{1}\mathrm{,}{k}_{2}}^{34},$ and ${\rho}_{{k}_{1}\mathrm{,}{k}_{2}}^{32}.$ Every loop containing this region, i.e. *f*_{23}, *f*_{13}, *f*_{24}, and *f*_{14}, must be expanded into combinations of diffraction orders both coming back and forth. In , an example of the different loop combination of *f*_{23} is given for orders –1, 0, +1.

In this calculation, a sufficient number of diffraction orders and guided modes must be taken into consideration. Whenever $\sqrt{{\u03f5}_{i}}{k}_{w}\mathrm{sin}\mathrm{(}\theta \mathrm{)}+\frac{m2\pi}{d}<{k}_{w}\sqrt{{\u03f5}_{i}},$ the diffraction order inside the medium becomes evanescent, and so the term *e*_{h} will become significantly small. As for the guided modes, as only the fundamental mode for the metallic slits couples strongly to the plane waves of regions 1, 3, and 5, only this mode is considered for the remainder of the article. The loop factors, with diffraction order (–1, 0, +1) considered for region 3 and only the transverse electric and magnetic (TEM)-guided mode considered in regions 2 and 4 (${\rho}_{\mathrm{|}0\mathrm{\u3009}\mathrm{,}|0\u3009}^{23}$ is written *ρ*^{23}), can be expressed as

$$\begin{array}{c}{f}^{12}={\rho}^{23}{\rho}^{21}{e}_{{h}_{1}}^{2}\text{\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}}{f}^{34}={\rho}^{45}{\rho}^{43}{e}_{{h}_{3}}^{2}\\ {f}_{{k}_{1}\mathrm{,}{k}_{2}}^{23}={\rho}_{{k}_{1}\mathrm{,}{k}_{2}}^{34}{\rho}_{{k}_{2}\mathrm{,}{k}_{1}}^{32}{e}_{{h}_{2}}\mathrm{(}{k}_{1}\mathrm{)}{e}_{{h}_{2}}\mathrm{(}{k}_{2}\mathrm{)}\\ {f}_{{k}_{1}\mathrm{,}{k}_{2}}^{13}={\rho}_{{k}_{1}\mathrm{,}{k}_{2}}^{34}{\tau}_{{k}_{2}}^{32}{\rho}^{21}{\tau}_{{k}_{1}}^{23}{e}_{{h}_{1}}^{2}{e}_{{h}_{2}}\mathrm{(}{k}_{2}\mathrm{)}{e}_{{h}_{2}}\mathrm{(}{k}_{1}\mathrm{)}\\ {f}_{{k}_{1}\mathrm{,}{k}_{2}}^{24}={\rho}^{45}{\tau}_{{k}_{1}}^{43}{\rho}_{{k}_{1}\mathrm{,}{k}_{2}}^{32}{\tau}_{{k}_{2}}^{34}{e}_{{h}_{2}}\mathrm{(}{k}_{2}\mathrm{)}{e}_{{h}_{2}}\mathrm{(}{k}_{1}\mathrm{)}{e}_{{h}_{3}}^{2}\\ {f}_{{k}_{1}\mathrm{,}{k}_{2}}^{14}={\rho}^{45}{\tau}_{{k}_{1}}^{43}{\tau}_{{k}_{1}}^{32}{\rho}^{21}{\tau}_{{k}_{2}}^{23}{\tau}_{{k}_{2}}^{34}{e}_{{h}_{2}}\mathrm{(}{k}_{2}\mathrm{)}{e}_{{h}_{2}}\mathrm{(}{k}_{1}\mathrm{)}{e}_{{h}_{1}}^{2}{e}_{{h}_{3}}^{2},\end{array}$$(8)where *k*_{1} is the parallel wavevector of light going forward and *k*_{2} is the one going backward, ${e}_{{h}_{2}}\mathrm{(}k\mathrm{)}={e}^{i{q}_{z}{h}_{2}}$ represents the propagation term inside the dielectric layer, and ${e}_{{h}_{j}}={e}^{i{k}_{{k}_{{z}_{j}}}{h}_{j}}$ the one for the metallic layers. The direct transmission of light in the structure is represented by

$${t}_{k}^{d}={\tau}^{12}{\tau}_{k}^{23}{\tau}_{k}^{34}{\tau}^{45}{e}_{{h}_{1}}{e}_{{h}_{2}}\mathrm{(}k\mathrm{)}{e}_{{h}_{3}},$$(9)with *k* wavevector inside the dielectric layer. The total transmission coefficient is given by

$$t={\displaystyle \sum _{k}{t}_{k}^{d}\mathrm{(}1+{A}_{14}\mathrm{)}},$$(10)where *A*_{14} represents the considered summation of every possible combination of loops of the system. This is achieved using the following matrix form:

$${A}_{14}={V}_{k}^{in}\times {[I-{T}_{k}]}^{-1}\times {V}_{k}^{out},$$(11)where ${V}_{k}^{in}$ is a vector representing the initial loop of a certain combination, **I** is the identity matrix, **T**_{k} is a matrix of a subsequent loop with the corresponding connecting coefficient, and ${V}_{k}^{out}$ represents the exit coefficient of the final loop considered. The middle term in Eq. (11) is a geometric series of matrices that sums up all possible combination of paths, in the same way a simple geometric series sums up all paths of a single layer. These vectors and matrices take the form

$${V}_{k}^{in}=\mathrm{[}{a}^{12}{f}^{12}{a}_{k\mathrm{1,}k2}^{23}{f}_{k\mathrm{1,}k2}^{23},\dots ,{a}_{k\mathrm{1,}k2}^{14}{f}_{k\mathrm{1,}k2}^{14t}\mathrm{]},$$(12)$${T}_{k}=\left[\begin{array}{cccc}{a}^{12\text{\hspace{0.17em}}-\text{\hspace{0.17em}}12}{f}^{12}& {a}_{k\mathrm{1,}k2}^{12\text{\hspace{0.17em}}-\text{\hspace{0.17em}}23}{f}_{k\mathrm{1,}k2}^{23}& \dots & {a}_{k\mathrm{1,}k2}^{12\text{\hspace{0.17em}}-\text{\hspace{0.17em}}14}{f}_{k\mathrm{1,}k2}^{14t}\\ {a}^{23\text{\hspace{0.17em}}-\text{\hspace{0.17em}}12}{f}^{12}& {a}_{k\mathrm{1,}k2}^{23\text{\hspace{0.17em}}-\text{\hspace{0.17em}}23}{f}_{k\mathrm{1,}k2}^{23}& \dots & {a}_{k\mathrm{1,}k2}^{23\text{\hspace{0.17em}}-\text{\hspace{0.17em}}14}{f}_{k\mathrm{1,}k2}^{14t}\\ \vdots & \vdots & \ddots & \vdots \\ {a}^{14\text{\hspace{0.17em}}-\text{\hspace{0.17em}}12}{f}^{12}& {a}_{k\mathrm{1,}k2}^{14\text{\hspace{0.17em}}-\text{\hspace{0.17em}}23}{f}_{k\mathrm{1,}k2}^{23}& \dots & {a}_{k\mathrm{1,}k2}^{14\text{\hspace{0.17em}}-\text{\hspace{0.17em}}14}{f}_{k\mathrm{1,}k2}^{14t}\end{array}\right],$$(13)$${V}_{k}^{out}={[\begin{array}{cccc}{t}_{{k}_{2}}^{12}& {t}_{{k}_{2}}^{23}& \dots & {t}_{{k}_{2}}^{14}\end{array}]}^{T},$$(14)where ${a}_{k\mathrm{1,}k2}^{xy\text{\hspace{0.17em}}-\text{\hspace{0.17em}}ab}$ represents a connecting coefficient that uses the previous wavevector *k*_{1} as the incoming wavevector for the reflection coefficient of the next loop. The exit coefficient ${t}_{{k}_{2}}^{23}$ inputs the correct transmission coefficient for the considered path, using the wavevector of the last loop. It is to be noted that the vectors and matrices of Eqs. (12) and (13) use the loop factors *f*_{13t} , *f*_{24t} , and *f*_{25t} , which consider the loop with a summation of the previous possible loop that can exist inside. Finer details are provided in the Supplementary Information section. The total transmission of the system is given by *T*=|*t*|^{2}. A Matlab code for the implementation of the model is available in the Supplementary Information section.