Multimodal Solution for a Circular Waveguide Radiating Into Multilayered Structures Using the Axially Symmetric Modes

In this article, an exact formulation is derived for the mode reflection coefficient (and in general, the full S-parameter matrix) for a circular waveguide radiating into a multilayered structure and when excited with any combination of axially symmetric modes (i.e., modes of the form TE0m and TM0m, where $\mathbf {m}$ is a positive integer). This formulation solves for the fields in the waveguide, including fields resulting from higher-order evanescent modes, using Fourier analysis. This leads to an accurate calculation of the mode S-parameter matrix, which includes the reflection coefficient of each excited mode. The derived formulations were validated through comparison to those computed using full-wave 3-D electromagnetic simulations. Additional simulations demonstrated the effect of considering higher-order modes on the results. The effect of having a finite flange and sample size on the complex reflection coefficient was also shown through 3-D simulations, indicating the relative insensitivity of the TE01 probe to edge effects. Reflection coefficient measurements, using a circular waveguide probe with the TE01 mode, were performed to experimentally verify the accuracy of the formulations. Forward-iterative optimization (i.e., optimal curve fitting) techniques were then performed on the reflection coefficient measurements to demonstrate the efficacy of this method for accurately estimating the thickness and complex permittivity of thin dielectric layers.


I. INTRODUCTION
O PEN-ENDED waveguide probes are widely used in many microwave and millimeter-wave nondestructive testing (NDT) applications, ranging from dielectric property measurement [1], [2], [3], [4], [5], thickness measurement [6], [7], porosity estimation [8], and detection of surface breaking cracks in metals [9], among others. Openended waveguide methods are particularly useful for the inspection of multilayered composite structures, given that information about both physical and electrical properties of each layer can be obtained [10], [11]. Using this method, the complex reflection coefficient, over a band of frequencies, is measured referenced to the waveguide aperture while radiating into a multilayered structure of interest [10], [11]. The measured complex reflection coefficient can then be used to solve for one or more parameters of the structure by using a forward-iterative optimization process along with a robust forward electromagnetic (EM) model that calculates the same reflection coefficient [10], [11].
In this approach, the choice of the forward EM model used in the forward-iterative procedure is critically important. While full-wave 3-D numerical EM solvers can be used as the forward model, they commonly suffer from issues related to both accuracy and long computation time. Since determining the structure parameters is done using an iterative optimization procedure, the forward model must be executed a large number of times, and thus the long computation time associated with numerical full-wave EM solvers becomes prohibitive. Many models have been developed that addressed this issue by providing approximate formulations that are quickly computed [12], [13], [14]. However, these formulations do not consider important factors, such as the influence of the higher-order evanescent modes that are present at the waveguide aperture, on the resulting complex reflection coefficient. As reported in [10] and [11], the error introduced by not considering these higher-order modes can be significant when calculating the desired parameters (electrical or geometrical) of the structure. To address this issue, an exact and rigorous formulation for the reflection coefficient measured by an open-ended rectangular waveguide excited with the dominant TE 10 mode and radiating into a generally multilayered structure was developed in [11], which considers the presence of the generated higher-order modes. This formulation is full wave and converges to an exact solution while still remaining computationally efficient.
While open-ended rectangular waveguide probes have been widely studied for the purpose of evaluating properties of multilayered structures, open-ended circular waveguides have not been studied nearly to the same extent. This is in part due to the fact that when operating in the dominant TE 11 mode, circular waveguides provide for less bandwidth than the rectangular waveguide and provide no other advantages. However, circular waveguides can provide unique advantages when operating in modes other than the dominant TE 11 mode. For example, the axially symmetric modes (i.e., modes of the form TE 0m and TM 0m for integer m) have been utilized in several applications, such as conductivity and loss factor measurement [15] and polarimetric SAR imaging [16], [17]. As will be shown in this article, operating an open-ended circular waveguide using these modes can also be useful for multilayered structure characterization and can provide some distinct and practically useful advantages over their rectangular waveguide counterparts. Specifically, this method is much more accurate at estimating the thickness and dielectric properties of thin and conductor-backed dielectric structures. This improvement in accuracy is largely due to the insensitivity of this method to sample and flange edges when using the TE 0m modes, which is a known problem with the rectangular waveguide method [11].
In this article, a rigorous and exact formulation is derived for the reflection coefficient measured as seen by an openended circular waveguide radiating into a multilayered structure while it is excited with any combination of axially symmetric modes. This formulation was derived using a similar approach to that used in [11] and accounts for the effect of higher-order modes generated at the aperture. Use of this formulation provides for the same accuracy as full-wave 3-D simulation tools but is much more computationally efficient (typically more than 10 000 times faster), as will be shown in this article. Additionally, the formulation provides the full-mode S-parameter matrix (i.e., reflection coefficients for each mode), which can potentially be used alongside a multimodal probe which supports more than one propagating mode.

II. FULL-WAVE COMPLEX REFLECTION COEFFICIENT FORMULATION
In this section, an expression for the complex reflection coefficient of an open-ended circular waveguide, excited with any axially symmetric mode, radiating into a general multilayered structure is derived. The axially symmetric modes are the circular waveguide modes of the form TE 0m and TM 0m , where m is a positive integer. This formulation is full wave and accounts for both propagating and evanescent modes and allows for computing the mode scattering (S-parameter) matrix. Fig. 1 shows the geometry of the problem. An open-ended circular waveguide, with an infinitely extended flange made of a perfect electric conductor (PEC) radiates into a multilayered dielectric structure. The layered structure consists of any number of stacked uniform infinite sheets (i.e., layers), each of which is described by its thickness, complex permittivity, and complex permeability. The last layer is either an infinite half space of a material or a PEC layer (i.e., conducting substrate). Note that the geometry of the problem is axially symmetric about the z-axis. Consequently, as long as the fields are excited with the same symmetry (e.g., TE 0m and TM 0m modes), the problem reduces to a 2-D one (i.e., z and ρ are the only relevant dimensions).
Note that the multilayered medium in Fig. 1 is described using relative complex permittivity ( r = / 0 ) and permeability (μ r = μ/μ 0 ). Complex permittivity can be separated into its real and imaginary parts (i.e., r = r −j r ), where r describes the ability of the medium to store electrical energy, and r is the loss factor of the medium. Additionally, the ratio r / r is referred to as the dielectric loss tangent of the medium.
In the past, a rigorous and exact formulation for the reflection coefficient of an open-ended rectangular waveguide radiating into a multilayered structure has been derived [11]. Here, we perform similar derivations for an open-ended circular waveguide excited with the axially symmetric modes. While the overall procedure for doing so is similar to that in [11], there are sufficient and important differences which warrant showing details of the derivation. For example, the fact that the geometry of the entire problem has axial symmetry (i.e., no dependence on φ) reduces this to a 2-D problem, as opposed to the 3-D problem considered in [11].
In this article, the magnetic and electric Hertzian vectors ( h and e , respectively) are utilized in order to simplify the derivations. The relationship between the Hertzian vectors and the fields are given explicitly in (1) and (2), where k 0 is the free-space wavenumber and η 0 is the free-space impedance One important fact to note is that the axially symmetric nature of this problem means that the TE 0m and TM 0m modes are completely independent of each other (i.e., there is no coupling between the TE 0m and TM 0m modes). Consequently, regardless of the specific parameters of the multilayered structure, any excited TE 0m mode only generates other TE 0m modes (i.e., no TM 0m modes are generated). The same can be said for the TM 0m modes. In addition to simplifying the derivations significantly, this observation allows for the formulation for the TE 0m and TM 0m modes to be completely separated. Thus, in this article, two formulations are considered and outlined, one for the TE 0m modes, and one for the TM 0m modes.

A. SOLUTION FOR TE 0M MODES
Before deriving the formulations in this section, several key variables and functions must be defined. Specifically, k c,0n = p 0m /a is the cutoff wavenumber of the TE 0m mode, J v (x) is the vth-order Bessel function of the first kind, p 0m is the mth root of J 1 (x), k 0n = k 2 0 r wg μ r wg − k c,0n 2 is the propagation constant of the TE 0m mode, k 0 = (2π/λ) is the free-space wavenumber, a is the radius of the circular waveguide, and r wg and μ r wg are the complex relative permittivity and permeability of the material filling the waveguide (if any). Let us consider that the waveguide in Fig. 1 is excited with a combination of TE 0m modes. The excited (or incident) fields can then be written as a superposition of the field distribution for each of these modes. Equation (3) gives the z-component of the magnetic Hertzian vector describing the superposition of these modes (the ρand φ-components are both zero). Deriving the E-field and H-field vectors, using (1) and (2), and confirming that they correspond to a superposition of TE 0m modes with mode coefficients A i,h 0m is fairly straightforward. We can similarly define the reflected modes using the reflected mode coefficients A r,h 0m , as done in (4). As mentioned earlier, the axial symmetry of this problem allows for neglecting all other modes, as there is no mechanism for generating them. Thus, the total fields in the waveguide are represented exactly by i,h z (ρ, z)+ r,h z (ρ, z). Note that, in (3) and (4), the mode coefficients are scaled such that if A i,h 0m = A i,h 0n , the transmitted power level for the TE 0m and TE 0n modes is the same The goal of these derivations is to provide a formulation for the mode S-parameter matrix referenced to the aperture of the waveguide. The S-parameters can be found by taking the ratio is zero for all q = n, and S mn is the S-parameter corresponding to transmitting the TE 0n modes and receiving the TE 0m mode. Thus, if a formulation for the reflected mode coefficients A r,h 0m given any arbitrary set of incident mode coefficients A i,h 0m can be found, then the computation of S-parameters becomes quite straightforward.
The fields in the multilayered structure can be written as a superposition of the forward and backward traveling waves (i.e., the waves traveling in the +zand −z-directions, respectively). Equation (5) gives the z-component of the magnetic Hertzian vector in the first layer expressed as a zeroth-order Hankel transform of the unknown spectral functions A h+ (k ρ ) and A h− (k ρ ), where the definition of a Hankel transform is given in (6). It should be noted that due to the symmetry in this problem, only the z-component of the magnetic Hertzian vector can be nonzero given the TE 0m mode excitation. Note that we are utilizing the fact that the zeroth-order Hankel transform of a function is equivalent (after scaling) to the 2-D Fourier transform of an axially symmetric function Additionally, it is worth noting that the spectral functions A h+ (k ρ ) and A h− (k ρ ) are the same as those defined in [11]. In this article, an analytical formulation for the reflection coefficient looking into the first layer ( h 1 (k ρ ) = A h− (k ρ )/A h+ (k ρ )) is derived, which we make use of here.
Note that k z 1 = k 2 0 r 1 μ r 1 − k 2 ρ , and that k ρ is the spectral coordinate corresponding to the spatial coordinate ρ. As mentioned in [11], it is relatively straightforward to verify that the fields defined by (5) satisfy the wave equation.
Next, we can find the relationship between the tangential electric and magnetic fields spectra in the first layer at the aperture (i.e., at z = 0). Combining (1) and (2) with (5) and applying the definition of the Hankel transform gives the following expressions for E t φ and H t ρ in the first layer (note that E ρ and H φ are zero everywhere): Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply.
Using (7) and (8), the tangential field spectra at the aperture can be found as shown in (9) and (10). Note that here, the spectra of E φ and H ρ are derived using the first-order Hankel transform, as opposed to the zeroth order. While the first-order Hankel transform is not equivalent to the Fourier transform, it is still invertible. Thus, as long as the first-order transform is consistently used for all field spectra, there is no issue. In this article, the hat accent is used to describe the spectrum of a spatial function, i.e., Evaluating the Hankel transforms in (9) and (10), noting that the Hankel transform is its own inverse, and combining and rearranging gives the relationships in (11) and (12). Note that (12) is expressed in terms of the TE reflection coefficient looking into the first layer, i.e., h , which as mentioned earlier is known analytically [11] We can similarly solve for the electric field spectrum inside the waveguide at the aperture. Similarly, we define the electric field spectrum to be the first-order Hankel transform of the electric field, as seen in (13) and (14). Note that E i φ (ρ, 0) and E r φ (ρ, 0) can be found by combining (1) with (3) and (4), respectively, and evaluating at z = 0. Also, note that the electric field is zero for ρ > a due to the PEC flange Applying the definition of the Hankel transform gives a direct formula for E wg φ (k ρ , 0) in terms of the mode coefficients, as seen in Next, we can enforce the magnetic field boundary conditions at the aperture, as done in (17). Taking the Hankel transform of both sides gives the relationship in (18) We can then multiply both sides of (18) by ρJ 1 (ρk c,0m ) and integrate over the aperture of the waveguide, as shown in The above expression can be simplified by evaluating the integral over the waveguide aperture and applying the definition of the Hankel transform. Note that since integration is a linear operation, the order of the Hankel transforms and the integral over the aperture can be interchanged. Using this fact and after some simplifications, the relationship in (20) is found The final step here is to apply the relationship between H t ρ (k ρ , 0) and E t φ (k ρ , 0) expressed in (11), and then enforce the electric field boundary conditions ( E t φ (k ρ , 0) = E wg φ (k ρ , 0)). After some simplifications, canceling of common terms, etc., this results in the relationship in The above expression can be used to calculate S-parameters, as will be described later. Note that further simplification of (21) is performed when showing the final equations in Section II-C.

B. SOLUTION FOR TM 0M MODES
Similar to what was done in the previous section, several variables and functions used throughout are defined. Specifically, k c,0n = p 0m /a is the cutoff wavenumber of the TM 0m mode, p 0m is the mth zero of J 0 (x), and k 0n = k 2 0 r wg μ r wg − k c,0n 2 is the propagation constant of the TM 0m mode. In this section, the derivation from Section II-A is repeated for TM 0m modes. Since the process is almost exactly the same, the full detailed derivation will not be shown. Instead, an overview will be given with details for only a few critical steps. The fields in the waveguide can now be expressed as a combination of TM 0m modes. As a result, the electric Hertzian potentials are used instead of the magnetic. Equations (22) and (23) give the z-components of the electric Hertzian vectors, where A i,e 0m and A r,e 0m are the excitation coefficients for the incident and reflected TM 0m modes. Again, due to the symmetry of the problem, no other modes need to be considered, and thus the total fields in the waveguide can be derived from the total electric Hertzian potential The electric Hertzian potential describing, the fields in the multilayer structure ( t,e z (ρ, z)), can be expressed in exactly the same way as in Section II-A, except that the unknown spectral functions areÂ e+ (k ρ ) andÂ e− (k ρ ), as shown in The relationship between the tangential electric and magnetic field spectra can be found using the same procedure as in Section II-A, resulting in (25) and (26). This relationship is expressed in terms of the TM reflection coefficient looking into layer 1 (i.e., h Finally, we can enforce the electric and magnetic field boundary conditions and perform the integration over the waveguide aperture, resulting in the relationship in Note that further simplification of (27) and (28) is performed when showing the final equations in Section II-C.

C. FINAL FORMULATIONS AND DISCUSSION
Based on the derivations in Sections II-A and II-B, the following formulations can be arrived at. For TE 0m mode excitation, the full formulation is given in k c,0n k c,0m I h mn (29) Similarly, the full formulation for TM 0m mode excitation is given in In the above formulations, the spectrum functions h 0 (k ρ ) and e 0 (k ρ ) were earlier defined in (12) and (26), respectively, which in turn are functions of the TE and TM reflection coefficient spectra looking into layer 1 [i.e., h 1 (k ρ ) and e 1 (k ρ )]. As mentioned earlier, these functions are given analytically in [11] and are repeated here in The formulations in (29)-(32) allow for the computation of the reflected fields at the aperture of the circular waveguide, given any axially symmetric excitation field definition, in turn allowing S-parameters to be computed. This can be performed by assembling and numerically solving the system of linear equations with respect to the mode coefficients. This is made easier by first expressing (29)-(32) in the form of a matrix equation containing the S-parameter matrix. This derivation is given in the Appendix-A and makes the implementation of these formulations much more straightforward. However, there are several caveats to consider. First, note that the formulation expresses the fields at the aperture using an infinite summation of mode coefficients. While this is not an issue when defining the incident modes, as typically only one mode is excited, there are in general an infinite number of reflected modes with nonzero magnitudes, even when there are only a finite number of propagating modes. Since we cannot numerically evaluate an infinite summation, we must use a truncated system of equations, only considering a finite number of modes. This technique works because the summations eventually converge to some specified accuracy. However, the question of how many modes to consider becomes important. This is explored in the next section.
The integrals in (30) and (32) cannot be evaluated analytically, and instead must be computed using numerical methods. Additionally, when evaluating lossless multilayered structures, the integrand will have singularities at several points. Similarly, when evaluating low-loss structures, the integrand will be close to singular at some points. To address this, adaptive numerical integration methods, such as global adaptive quadrature, must be used to compute the integrals within some convergence tolerance. However, compared to the formulations for rectangular waveguides [11], the process is significantly simpler and faster, as the integrals are 1-D instead of 2-D.
Because the S-parameter formulations given in this article are not analytically invertible with respect to the structure parameters, numerical forward-iterative optimization techniques must be applied (i.e., curve fitting). This optimization issue is discussed in detail in [10] and [11], and for the most part, the findings directly apply here. However, the use of circular waveguide axially symmetric modes provides a couple of advantages over the use of rectangular waveguide modes. First, the integrals that must be evaluated here are 1-D (i.e., integration with respect to k ρ only), as opposed to the two-dimension integrals required for the rectangular waveguide. The corresponding increase in computational efficiency means that the forward model can be evaluated significantly faster (more than 100 times faster), and thus the use of optimizers, especially global optimizers, is not so time consuming as to be prohibitive. Second, the use of multimodal S-parameter measurements (using a multimodal probe) allows for additional measurement diversity in information about the multilayered structure. This may allow for higher accuracy when solving for structure parameters. For example, it is possible to measure a structure using both a TE 01 mode probe as well as a TM 01 mode probe. The forward-iterative procedure can then be applied to both sets of measurements simultaneously, potentially providing additional information in the same way that frequency diversity is currently used.

A. EFFECT OF HIGHER-ORDER MODES
In order to validate the formulations derived in Section II, the calculated complex reflection coefficient seen by a circular waveguide radiating into a layered structure and excited with either the TE 01 mode or the TM 01 mode was compared to simulations performed by the full-wave 3-D numerical EM solver available in CST Microwave Studio.
As mentioned earlier, the formulations in this article express the reflection coefficient of interest using an infinite summation of mode coefficients. In order to calculate these summations numerically, the summations must be truncated such that only a finite number of modes are considered. Even though in this case only the TE 01 mode is propagating, consideration of the presence of nonpropagating (i.e., evanescent) modes is important. To demonstrate this fact, Fig. 2(a) compares the frequency-dependent complex reflection coefficient of the incident TE 01 circular waveguide radiating into the air as calculated using CST Microwave Studio (designated as 3-D Simulations in Fig. 2) and using the formulations in this article. As expected, as the number of considered modes increases, the results approach that  of the full-wave EM solver (i.e., CST Microwave Studio). Fig. 2(b) shows a similar plot for the TM 01 mode, where the error is slightly larger than for the TE 01 mode. Note that the frequency range used in these plots is 1.25f c -1.75f c , where f c is the cutoff frequency of either mode. Although, based on the results in Fig. 2, one might conclude that only a specific number of modes may need to be considered to obtain a desired accuracy, the specific number of required modes depends on a variety of factors, namely, the waveguide dimensions, the operating frequency, and the properties of the multilayered structure. To demonstrate this important fact, Fig. 3 illustrates the same comparison of complex reflection coefficient but for a waveguide looking into a conductor-backed layer with r = 4−j0.05 and a thickness of λ c /32, where λ c is the cutoff wavelength of either of the propagating modes. As can be seen, the TE 01 mode reflection coefficient remains almost unchanged regardless of the number of modes considered, while the TM 01 mode reflection coefficient is highly dependent on the number of modes (due to the increased magnitude of the higher-order modes), more so than in Fig. 2. It should also be noted that the accuracy and convergence settings of the 3-D numerical EM solver had to be set extremely high to get results that matched the formulations in this article.
A similar test was performed to verify the accuracy of the formulation for multilayered structures. Specifically, a 3-layer structure was considered consisting of a λ c /8-thick layer of air representing a standoff, a λ c /2-thick layer with r = 4 − j0.05, and finally an infinite half-space layer of  air. Fig. 4 shows the calculated reflection coefficients for either mode when considering various number of modes and compared to the results found using the 3-D full-wave solver. Clearly, as the number of considered modes increases, the calculated reflection coefficient approaches that of the 3-D numerical EM solver.
The use of the formulations in this article allows the calculation of reflection coefficient with much higher computational efficiency for the same level of accuracy. As an example, calculation of the TE 01 mode reflection coefficient in Fig. 4(a) at 51 uniformly spaced frequency points and with an accuracy of −60 dB required ∼1425 s when using CST Microwave Studio, and ∼0.06 s when using the formulations in this article (with 15 modes considered), which is a speedup of ∼23 750 times. Similarly, the single-layer structure in Fig. 3(a) required ∼88 s when using CST Microwave Studio, and ∼0.005 s when using the formulations in this article (with three modes considered), resulting in a speedup of ∼17 600 times.
Another issue that should be considered is the presence of resonant higher-order modes. When a TE 01 mode radiates into a low-loss (e.g., loss tangent less than 0.0001), high permittivity/permeability (e.g., r > 20) conductorbacked one-layer structure having an electrical thickness less than but approaching λ/2, resonances can occur. This effect was previously reported for the TE 01 mode, and in fact, the behavior of the reflection coefficient at the resonance frequency was used to characterize conductivity and loss factor with high accuracy [15].
To illustrate this effect, Fig. 5(a) and (b) shows the TE 01 mode reflection coefficient versus frequency for a conductorbacked dielectric layer structure with r = 25 − j0.001 and a thickness of 0.045λ c . Note that in this simulation, five modes are considered; otherwise, these resonances do not all manifest themselves. These are caused by a sharp increase in the magnitude in each of the evanescent modes, which can be confirmed by Fig. 5(c), which shows the magnitude of each of the first five mode coefficients versus frequency. Note that if only one mode is considered, the reflection coefficient is exactly the same except that the resonances are no longer present. Although it may seem that the presence of these resonances could be caused by an error in the formulations , the correctness can be confirmed using a full-wave numerical solver (CST Microwave Studio). To demonstrate, Fig. 6(a) shows the electric field magnitude for this case inside the waveguide and the dielectric layer at 1.564f c , which corresponds to the third from the left resonance seen in Fig. 5(b) and (c) and the second from the top resonance seen in Fig. 5(a). The electric field distribution at the aperture is clearly dominated by the TE 04 mode, as expected from Fig. 5(c). In contrast, Fig. 6(b) shows the same distribution but at 1.55f c , which is slightly off the resonance. In this plot, the TE 01 mode dominates, which is the normal case.
From these simulations, we can conclude that the number of considered modes must be changed on an applicationby-application basis, as there are several different factors involved by which to give a general rule. Alternatively, the number of modes can be adaptively changed for each calculation by increasing the number of modes until the desired convergence criteria are reached.

B. EFFECT OF SAMPLE/FLANGE EDGES
The formulations derived in this article assume that the waveguide flange and the multilayer structure both extend infinitely in the x-and y-dimensions. However, this assumption is not true in practice, as most physical geometries of interest are finite in extent (in some way). For rectangular VOLUME 2, 2023 2500208  waveguides, the effect of a finite flange dimensions has been reported to be significant, especially for conductor-backed, thin structures [11]. In fact, significant effort has been expended to remedy this issue by using specially designed flanges [18]. In this section, the effect of the finite flange and the sample edge is investigated. To this end, CST Microwave Studio was used to simulate the reflection coefficient as a function of flange/sample size for a one-layer conductor-backed structure. Fig. 7 shows the relevant geometry, where the flange is infinite in all directions except for the +y-direction. The distance between the edge of the flange and the closest edge of the circular waveguide is denoted by the distance d. The sample is also truncated in the same manner as the flange (i.e., the flange edge and the sample edge have the same x-and y-coordinates).  Fig. 8 shows the reflection coefficient error due to finite dimensions (or finite flange error), which is calculated as err = | inf − finite |, where inf is the reflection coefficient calculated assuming an infinite flange and sample, and finite is the simulated reflection coefficient when the flange and sample are finite. This error is plotted for both the TE 01 and TM 01 modes and for two different structures. Specifically, both structures have a single layer of dielectric with r = 4 − j0 but the structure thickness is either λ/8 or 3λ/8 and when conductor-backed. It is clear that the error is higher for both modes for the 3λ/8 layer. Note that for the TE 01 mode, the error is nearly zero for the λ/8 layer, even when the edge of the flange is almost touching the waveguide.
To obtain some insight into why the finite flange error is small when using the TE 01 mode with thin structures, simulations were performed showing the electric field distribution in this case. Fig. 9(a) shows the magnitude of electric field distribution inside a circular waveguide exited with the TE 01 mode and radiating into a single-layer conductor-backed dielectric with r = 4 − j0 and a thickness of λ/8. Immediately we can notice that the fields within the onelayer structure are contained within a region close to the waveguide aperture. The reason for this can be found by noting that the flange conductor and the backing conductor effectively form a parallel plate waveguide. In order for the waves to interact with the flange and sample edges, they must propagate within the parallel plate waveguide. However, when the waveguide is excited with the TE 01 mode, only TE modes are present (i.e., the z-component of the electric field is always zero). As a result, any mode generated within the parallel plate waveguide is evanescent if the electrical thickness of the dielectric layer is less than λ/2 (since the cutoff frequency of any such mode will be higher than the operating frequency). Note that the layer electric thickness refers to the thickness relative to the wavelength in the dielectric (i.e., t √ r 1 μ r 1 ). This phenomenon explains why the finite flange error is larger for the 3λ/8 structure, as the electrical thickness is greater than λ/2. Additionally, the TM 01 mode does not have this property, as there exist propagating TM modes in the parallel plate waveguide regardless of the frequency, as shown in Fig. 9(b). Additional intuition about the difference between these two modes can be found by noting that the electric field magnitude vanishes at the wall of the waveguide for the TE modes, but not for the TM modes, leading to differences in the amount of coupling into the layered structure.
This feature of the TE 01 mode circular waveguide is beneficial in many practical situations. Specifically, it means that this method can be used on small samples or with a small flange, as will be demonstrated through measurement later. Additionally, the fact that the fields are contained within the region nearby the waveguide aperture means that the magnitude of the reflection coefficient directly corresponds to the absorption of the signal by one or more layers (i.e., due to a lossy medium), allowing for more accurate measurement of loss factor.

IV. MEASUREMENT RESULTS
In order to demonstrate the use of the formulations derived in this article, reflection coefficient measurements were performed using a circular waveguide excited with the TE 01 mode. To perform these measurements, a TE 01 mode feed for a circular waveguide was first designed and fabricated. Although many TE 01 mode circular waveguide feeds can be found in [19], [20], and [21], the feed used in this article was based on the design in [19], which was chosen since it provides excellent mode purity. Fig. 10(a) shows a view of the interior of the waveguide feed, which converts a TE 10 rectangular waveguide mode into a TE 01 circular waveguide mode. The feed structure consists of a rectangular waveguide which is split equally in four ways, each of which feeds into a circular waveguide. Fig. 10(b) shows the full view of the probe. This exact probe was also used in [15]. The probe is designed to be fed using a standard Ka-band (26.5-40 GHz) rectangular waveguide (WR-28) input, while the corresponding circular waveguide radius is ∼5.8 mm. Although the waveguide frequency operating range covers the entire Ka-band of frequencies, the operating band of the circular waveguide only includes the 32-40 GHz portion of the Ka-band (26.5-40 GHz) rectangular waveguide. The probe shown in Fig. 10 was fabricated from aluminum. A picture of this probe is shown in Fig. 11.
In order to accurately measure the reflection coefficient, a proper calibration of the probe must be performed. This can be performed using the standard 1-term network analyzer calibration technique by measuring three known calibration loads. However, at the time of writing this article, we had no access to a circular waveguide matched load, which is typically required for this purpose and for calibration accuracy. Instead, two of the circular waveguide feeds shown in Fig. 11 were fabricated, and a 2-port network analyzer was used to perform a thru-reflect-line (TRL) calibration [22], using a short circular waveguide section as the "line" calibration standard.
For all measurements shown in this article, the following setup was used. An Anritsu MS4644A vector network analyzer (VNA) was configured for two-port WR-28 waveguide measurements at 401 uniformly spaced frequency points covering the 32-40 GHz frequency range. Two coaxialto-waveguide adaptors (2.92 mm coaxial to WR-28) were used to connect the network analyzer ports to the feed ports of each probe. The TRL calibration method described above was used to calibrate with reference to the waveguide apertures. One of the waveguides was then used to make reflection coefficient measurements over the entire frequency band.
To demonstrate the utility of this method for thickness and complex permittivity measurement, the reflection coefficient when radiating into two different conductor-backed dielectric layer structures was measured. First, a 1.06-mm thick quartz ( r ≈ 4−j0.0004) microscope slide was measured. Fig. 12(a) shows a picture of the circular waveguide and flange with the quartz slide placed such that it is centered over the waveguide. Fig. 12(b) shows the same but with the quartz slide situated so that one of its edges is aligned with the edge of the waveguide. For either case, a flat aluminum metal plate was placed on top of the quartz slide, acting as the conductor backing. Then, the reflection coefficient of the probe was measured over the range of frequencies specified earlier. Fig. 13 shows the reflection coefficients plotted for both cases. Note that the reflection coefficient is almost completely unchanged, even for the extreme case, where the sample edge and the circular waveguide edge are touching. This demonstrates the measured insensitivity to the sample edge when using the TE 01 mode for thin conductor-backed structures, as expected. Note that the curve in either case is not smooth as would be expected. This is likely due to a calibration error.
The measured reflection coefficient data was then used to retrieve the complex permittivity of the layers via forwarditerative optimization of the model developed in this article. Specifically, a local optimizer was used to optimize complex permittivity and thickness such that the calculated reflection coefficients best matched the measurements (i.e., to minimize RMS error). The outputs of this optimization procedure are the real and imaginary parts of the complex permittivity along with the thickness. Note that here we assume that the complex permittivity does not change over the relatively small bandwidth used, allowing us to use the frequency diversity advantageously. Fig. 14 shows a plot of the measured reflection coefficient along with the best-fitting calculated reflection coefficient.
Using the results of the optimization procedure, the quartz was estimated to be 1.061-mm thick with a complex permittivity of r ≈ 4.02−j0.002. Note that the estimated thickness of the quartz is very close to the measured thickness of 1.07 mm. Additionally, while the real part of the estimated quartz complex permittivity (i.e., dielectric constant) is in the expected 3.8-4.3 range reported in [23] and [24], the loss factor is larger by a factor of about 10. However, the loss factor in this case is very difficult to measure due to its low value. This is as expected, since resonant cavity methods are typically required to measure very low loss factors (e.g., [15]). As will be shown, this method is much more accurate at measuring higher loss factors. The same forward-iterative procedure was used to estimate the thickness and complex permittivity using the measurements made with the quartz edge touching the waveguide opening. Using these measurements, the quartz was estimated to be 1.059-mm thick with a complex permittivity of r ≈ 4.02 − j0.002, which is very close to the earlier case, demonstrating the insensitivity of this method to sample edge effects.

A. COMPARISON TO RECTANGULAR WAVEGUIDE MEASUREMENTS
In order to demonstrate the advantages of using the circular waveguide TE 01 mode for the characterization of multilayered structures, measurements were performed and compared to those made with an open-ended rectangular waveguide. To perform a fair comparison, the engineered open-ended rectangular waveguide flange designed in [18] was used. This flange was designed such that reflection coefficient measurements made with it are as close as possible to measurements made with a theoretically infinite flange, as outlined in [18]. A picture of this flange designed for operation in the Ka-band (26.5-40 GHz) is shown in Fig. 15.
The flange picture in Fig. 15 was attached to the same VNA used in the previous section using a coaxial-towaveguide adaptor (2.92 mm coaxial to WR-28). The VNA was then configured for reflection coefficient measurements at 401 uniformly spaced frequency points covering the 26.5-40 GHz frequency range. Measurements of commercially available waveguide calibration standards were then used to calibrate the VNA with reference to the waveguide aperture.
Next, both the circular waveguide (using the setup described in the previous section) and the rectangular waveguide were used to make measurements of FR-4 conductorbacked sheets with two different thicknesses of 0.84 and 3.20 mm. These physical dimensions correspond to electrical thicknesses of ∼ λ/4 and ∼ λ, respectively, at 40 GHz. The reflection coefficients for either probe when looking into these two structures were measured and subsequently used to estimate the complex permittivity and thickness using the forward-iterative optimization procedure described earlier.  Note that for the rectangular waveguide, the forward model derived in [11] was used. For each probe, these measurements were repeated for a total of ten times, with each subsequent measurement being made after moving the probe to a different location on the sheet. Note that since the sheets are finite in the x-and y-directions, the probe was moved to several locations near an edge and several locations that were far from any edge, allowing us to compare sensitivity to edge locations. Table 1 shows the average and standard deviation of these estimated structure parameters for the 0.84-mm thick FR-4 sheet over the course of the ten measurements. These measurements were compared to those made with the filled waveguide technique, in which the samples are cut to fit into a rectangular waveguide section whose 2-port S-parameters are then measured [25] and used to solve for complex permittivity. This method is destructive yet very accurate, and thus is used as a reference measurement. Additionally, thickness measurements using a caliper are used as the reference thickness measurement.
From these results, we can see that the circular waveguide measurements yield more accurate results than the rectangular waveguide. Although the standard deviations for the rectangular waveguide measurements are low, there is a significant error in the estimated structure parameters, especially with regard to thickness and the loss factor.
The results from the measurements of the 3.20-mm thick FR-4 sheet are shown in Table 2. Again, the average value and standard deviation are shown for each parameter over the course of the ten measurements. In this case, the results for the circular waveguide are less accurate than for the thin sample, as expected. Specifically, the real and imaginary parts of the permittivity were overestimated, although the thickness and loss factor estimates were still much more accurate compared to the rectangular waveguide estimates.
Overall, the parameter estimation using the circular waveguide excited with the TE 01 mode was shown to be more accurate than when using the rectangular waveguide, even without the use of a specially designed flange.

V. CONCLUSION
In this article, a full-wave formulation was derived for the mode S-parameter matrix of a circular waveguide looking into a multilayered structure when excited with any axially symmetric modes (i.e., modes of the form TE 0m and TM 0m , where m is a positive integer). This formulation is exact and considers higher-order evanescent modes present at the waveguide aperture. Simulations were performed to validate the derived formulations through comparison with results using a full-wave 3-D EM solver. The use of this formulation allows accurate and computationally efficient calculation of the desired S-parameters (computation time is typically more than 10 000 times faster compared to general-purpose full-wave solvers). Additional simulations demonstrated the effect of considering higher-order modes on the results. The effect of having a finite flange and sample size on the reflection coefficient seen by the waveguide was shown through simulations using the 3-D numerical EM solver. The results of these simulations show that when using the TE 01 mode to interrogate thin, conductor-backed structures, the sensitivity to the flange and sample edges is almost nonexistent. This feature was confirmed experimentally through reflection coefficient measurements made with a circular waveguide exited with the TE 01 mode. Forward-iterative optimization using the developed formulations (i.e., curve fitting) was then performed on the measurements, demonstrating the efficacy of this method for the estimation of thickness and complex permittivity of thin dielectric layers.

APPENDIX
The formulations given in (29)-(32) are readily expressed in matrix form. Specifically, the mode S-parameter matrix of interest can be found by solving the following matrix equation. The S-parameter matrix can be written into the following matrix equation, where S is the mode S-parameter matrix, and the elements of A, B, K A , and K B can be found through inspection of either (29) and (30) or (31) and (32), depending on whether TE or TM modes are being solved. Clearly, once A, B, K A , and K B are known, S can be trivially solved for

A. MATRIX FORM FOR TE 0M MODES
Equations (38)-(41) define the elements of the square matrices A, B, K A , and K B used in (37). The size of each matrix is determined by the number of modes considered. In this case, the element S mn corresponds to the reflection/transmission coefficient seen by the circular waveguide when transmitting with the TE0n mode and receiving with the TE 0m mode. Note that B, K A , and K B are diagonal matrices and thus only their diagonal elements are specified