# Analytic Solutions to Reflection-Transmission Problem of Interface in Anisotropic Ice Sheet

Authors: Bangbing Wang; Honkuan Wong
Journal: International Journal of Antennas and Propagation
 (2020)

---

## Abstract

The rheology and evolution of the polar ice sheet are deeply influenced by the anisotropy of ice crystals. Studying the anisotropy of ice crystals can help to well understand and predict the behavior of the polar ice sheet and then the sea level rising and global climate change. In this paper, firstly, we deduce the expression of eigenvalues and eigenvectors of anisotropic media, which are determined by permittivity tensor and geometry of media. Then, the analytic formulas of reflection and transmission coefficients are derived directly by matrix transformation. Some models with real ice parameters are tested, and they present some special features at the anisotropic interface. We also discuss the physical meanings of eigenvalues and eigenvectors and the geometry analyzing to polarimetric radar. This analytic solution reveals the functional relationship between the macroradar reflection and the microphysical properties of ice crystals, which provides a feasibility of ice fabric identification by polarimetric radar detection.


---

## Body

## 1. Introduction

The polar ice sheets play an important role in the global climate system. Their evolution and stability have great impact on climate change and sea level rising [1, 2]. The huge thick ice sheet [3] (more than 2,000 m on average) is composed of numerous microice crystals. A single ice crystal shows noticeable anisotropy in its response to stress, specified by its exclusive C-axis, which is orthogonal to its basal plane. However, the bulk of ice crystals may exhibit preferred crystal orientations (fabric) and variable shapes and sizes (texture) under the giant stress in the ice sheet [4, 5]. The preferred fabrics and textures of ice crystals significantly affect the flow and evolution of ice [6–8]. Therefore, the fabric alignment in the ice sheet is a key indicator for understanding the behavior of the ice sheet. Unfortunately, this information can only be obtained from very few deep ice cores in Antarctic and Greenland [9–12].Glen [13] presented an isotropic ice flow law, which successfully explained the early field observations and had been widely accepted. However, recent field measurements showed some discrepancies with the results predicted by Glen’s flow law, given the improvements from equipment detection precision and extent [14]. These discrepancies implied that the anisotropy of the ice crystal could not be ignored in understanding the ice sheet evolution. Considering the importance of ice anisotropy in the ice sheet evolution [9, 14], how to recognize the fabric and its distribution in the ice sheet and how to understand its anisotropic behavior are still the leading and challenging fields in polar research [15].RES (or radio-echo sounding) has been widely employed as the standard facility in ice sheet expeditions since the 1950s. It achieves remarkable success in the polar exploration due to its high efficiency and accuracy [16, 17]. Early works since the 1950s had reported the elliptical polarization in RES echoes in ice sheet survey [18]. Hargreaves [19] suggested birefringence as the cause of echo difference along different polarized orientations and proposed an anisotropic expression of the permittivity with tensor. It is till the late 1990s that the permittivity tensor of the ice crystal within the frequency band of radar waves is precisely measured in the laboratory [20, 21]. Some field surveys have been conducted to address the relation between the polarization of the EM wave and the ice fabric in the ice sheet [22–26]. Additionally, some numerical methods have been developed for modeling the propagation of electromagnetic waves in the ice sheet. However, only a few [27–29] incorporated the anisotropy of the ice. Polarimetric radar is considered as the most potential tool for revealing the evolution of the anisotropic ice sheet. The major obstacle of using it in polar research is the limited understanding of radar wave propagation in the anisotropic ice sheet rather than hardware. This limited understanding deeply bottlenecks the data processing and interpretation of the polarimetric radar.The ice sheets can be simplified as a horizontally homogeneous stratified model with multiple anisotropic layers (isochronous layers) [27–29], which possess variable fabrics and permittivity tensors changed with the paleoclimate. Therefore, while considering the propagation of the EM wave in the ice sheet, the solution for this problem can be categorized into two sequential steps: the reflection and transmission at the anisotropic interface and the propagation in stratified and anisotropic media. For stratified anisotropic media, Berreman’s 4 × 4 propagation matrix is the most popular method used in optic and other EM fields [30–35]. Ursin [36] presented a propagation matrix method commonly used for elastic and EM wave propagating in horizontally layered media. Sluijter et al. [37] used the polarized ray-tracing method for analyzing the inhomogeneous anisotropic media. However, in ground penetrating radar, the wave sources are often transient pulses, not the same as the continuous wave used in EM and optic fields. The echoes of radar can be considered as the convolution response of a ray series (or geometrical optics series) [38] with the excitation wavelet. Each ray series is the combination of multiple reflections and transmissions at interfaces along the ray path. So, the reflection and transmission at the interface is a key basis for understanding the EM propagation in stratified media. There are some approaches to determine the reflections and transmissions of the EM waves at the interface of anisotropic media [39–42]. In this paper, we present a matrix transformation method to derive the direct expressions of reflection and transmission at the interface. It will be useful for the propagation in stratified and anisotropic media.


## 2. Methodology

### 2.1. Simplification of the Electromagnetic Model

RES survey and deep ice core analysis have proved that, in small scale, the deep ice sheet can be simplified as horizontal stratified ice layers. The ice sheet is often composed of four main fabric types [6, 27].Figure1 demonstrates the Schmidt diagrams with the four fabrics in the ice sheet. C-axis is a direction orthogonal to the crystal basal plane and often acts as the unique indicator of a crystal. The gray points or filled areas in Figure 1 indicate the possible C-axis directions or the gathering shape of C-axis (fabric) in polycrystalline ice. The symbols x and y indicate the principal axes in the horizontal plane, and the z-axis is perpendicular to the horizontal plane. Random ice (Figure 1(a)) often appears in the shallow part of the ice sheet because the stress is faint, and the C-axis points any direction kept while snow fell on the surface. The perfect single pole (Figure 1(b)) and vertical girdle (Figure 1(d)) are the two extreme states of fabric with ice crystal deformation under the stress and shear. The elongated single pole (Figure 1(c)) is the transient state between the perfect single pole and the vertical girdle under shear. Mathematically, anisotropic crystals in the ice sheet can be described by the second-order permittivity tensor, which has a diagonal form along the principal axis direction:(1)ε⟶=ε1000ε2000ε3.Figure 1
Schmidt diagrams demonstrate typical fabrics in the ice sheet: (a) a random ice, (b) a perfect single pole, (c) an elongated single pole, and (d) a vertical girdle.
(a)(b)(c)(d)The fabric can be categorized based on the signatures of differentiation of permittivity along principal axes. The permittivity of a crystal along or perpendicular to itsC-axis is denoted as ε⊥ and ε//, respectively, which had precisely been measured in the laboratory [20, 21]. The four main fabric types can be categorized through the permittivity tensor: a random ice (Figure 1(a)) has permittivity elements of ε1=ε2=ε3; a perfect single-pole ice (Figure 1(b)) has permittivity elements of ε1=ε2=ε⊥ and ε3=ε//; a vertical girdle fabric (Figure 1(d)) can be achieved when the permittivity elements satisfy ε1=ε// and ε2=ε3=ε//+ε⊥/2; and the permittivity tensor of an elongated single pole (Figure 1(c)) has ε1<ε2≪ε3.In radar analysis, EM wave propagations are typically described in two coordinates: measurement and media coordinates. Measurement coordinate is used to describe the amplitude and propagation direction of EM components. Media coordinate is used to define the principal axes of permittivity tensor and geometry of media. Figure2 presents the relative position of the measurement and media coordinates, where x′ denotes the semimajor axis of the indicatrix ellipses and y′ denotes the semiminor axis. Note that the z-axis in the measurement coordinate coincides with the z’-axis in the media coordinate because EM waves propagate as plane waves with normal incident angle into the ice sheet, where isochrones are generally horizontal. In the layered ice model, the measurement coordinate is shared among all layers, whereas the media coordinate differs in each layer and is simply a rotation of the measurement coordinate with an azimuth angle φaroundz/z’axis.Figure 2
Schematic diagram demonstrates the measurement coordinate ({x, y, z}) and media coordinate ({x′,y′,z′}). They own a common z/z' axis and rotate with azimuth angle φ. The red ellipse is the indicatrix ellipse, which indicates the principal axis direction of the permittivity tensor. Here, we suppose that the semimajor axis denotes x′.Applying the coordinate transformation to the diagonal permittivity tensor (formula (1)) yields the generic form of the permittivity tensor in the measurement coordinate:(2)ε¯=ε1sin2φ+ε2cos2φ12ε1−ε2sin2φ012ε1−ε2sin2φε1cos2φ+ε2sin2φ000ε3.


### 2.2. Eigenanalysis of EM Wave Propagations in Anisotropic Ice

The EM wave propagation in anisotropic ice-sheet layers can be described by differential-form Maxwell’s equations in time-harmonic field:(3)∇×E⇀i=−jωμ0H⇀i,∇×H⇀i=jωε¯iE⇀i,where E⇀i and H⇀i are the electric and magnetic fields in the i-th ice sheet layer, ω is the angular frequency, μ0 is the isotropic permeability of the i-th layer, and ε¯i is the permittivity tensor of the i-th layer. In the measurement coordinate, we suppose z-axis is the direction of EM propagation. It satisfies the case of normal incidence. So, electric field E and magnetic field H just exist in the plane of x- and y-axis determining. Combining electric and magnetic fields into one vector(4)ψi=ExEyHxHyTyields the matrix form of Maxwell’s equations [43–45]:(5)dψizdz=jω⋅Γiz⋅ψiz,where(6)Γi=000−μ000μ00ε1−ε2sinφcosφε1cos2φ+ε2sin2φ00−ε1sin2φ+ε2cos2φ−ε1−ε2sinφcosφ00.The entries ofΓi are expressed in terms of the azimuth angle φ and permittivity tensor in the i-th layer. For vertically propagating plane waves with wave number λ, the eigensolutions of (5) can be found by letting(7)ψ=ψ0eλz,which yields an eigenvalue problem with the eigenvalues λ:(8)Γi−λI⋅ψ0=0.In general, linear system (8) has four eigenvalues. Therefore, the general solution of (7) can be written as the linear combination of eigenvectors of (8):(9)ψz=A1a1ejλ1z+A2a2ejλ2z+A3a3e−jλ3z+A4a4e−jλ4z,or equivalently in matrix-vector notation:(10)ψz=a¯ejλ¯zA¯,where(11)a¯=a1a2a3a4,A¯=A1A2A3A4T,ejλ¯z=ejλ1zejλ2zejλ3zejλ4z.We can further reformulate (8) into(12)Γ−λI⋅X=0,where(13)λ⟶=λ1λ2λ3λ4=ε1μ0ε2μ0−ε1μ0−ε2μ0.Equation (12) has two fundamental solutions:(14)XI=x1x2x3x4=1ctgφλiμ0ctgφ−λiμ0,XII=x1x2x3x4=tanφ1λiμ0−λiμ0tanφ.For the normal incident EM wave along the -z/z′ axis, the wave field can be decomposed into waves I and II, where waves I and II are in the direction of x′ and y′ axis, respectively. In the measurement coordinate, waves I and wave II have rotation angles of φ and φ+π/2, respectively. The eigenvectors of the wave field in the i-th layer have two forms:(15)a¯i1=tgφi1tgφi11ctgπ2+φi1ctgπ2+φiεi1μ0εi2μ0⋅ctgπ2+φi−εi1μ0−εi2μ0⋅ctgπ2+φi−εi1μ0⋅tgφi−εi2μ0εi1μ0⋅tgφiεi2μ0,=tgφi1tgφi11−tgφi1−tgφi1ηi1−1ηi2⋅tgφi−1ηi11ηi2⋅tgφi−1ηi1⋅tgφi−1ηi21ηi1⋅tgφi1ηi2.(16)a¯i2=1tgπ2+φi1tgπ2+φictgφi1ctgφi1εi1μ0⋅ctgφiεi2μ0−εi1μ0⋅ctgφi−εi2μ0−εi1μ0−εi2μ0⋅tgπ2+φiεi1μ0εi2μ0⋅tgπ2+φi.=1−ctgφi1−ctgφictgφi1ctgφi11ηi1⋅ctgφi1ηi2−1ηi1⋅ctgφi−1ηi2−1ηi11ηi2⋅ctgφi1ηi1−1ηi2⋅ctgφi,where(17)ηi1=μ0εi1,ηi2=μ0εi2,are intrinsic impedances of waves I and II in the i-th layer, subscripts 1 and 2 represent waves I and II, respectively, and φi is the azimuthal rotation angle of the i-th layer in the measurement coordinate with respect to the media coordinate. Subscripts (1) and (2) indicate the two forms of matrices in (15) and (16), which are eigenvectors composed of tangent and cotangent functions of the azimuth angles, respectively. Taking the inverse of matrices (15) and (16) yields(18)a¯i1−1=12sec2φitgφi1tgφi11−tgφi1−tgφiηi1−ηi2⋅tgφi−ηi1ηi2⋅tgφi−ηi1⋅tgφi−ηi2ηi1⋅tgφiηi2T,(19)a¯i2−1=12csc2φi1−ctgφi1−ctgφictgφi1ctgφi1ηi1⋅ctgφiηi2−ηi1⋅ctgφi−ηi2−ηi1ηi2⋅ctgφiηi1−ηi2⋅ctgφiT.The inverse matrices (18) and (19) resemble the matrices (15) and (16), except for a transpose operator, inverse of wave impedance, and a scaling factor. Therefore, the inverse matrices (18) and (19) can be easily written once the dielectric tensor and azimuth angle are given.The calculation of reflection-transmission coefficients needs normalized eigenvectors, so we first normalize each column of the matrix in (15) by the L2 norm of the column vector:(20)Na¯i1=1secφitgφLi11Li2tgφiLi11Li21Li1−tgφiLi21Li1−tgφiLi21ηi1⋅Li1−tgφiηi2⋅Li2−1ηi1⋅Li1tgφiηi2⋅Li2−tgφηi1⋅Li1−1ηi2⋅Li2tgφiηi1⋅Li11ηi2⋅Li2,where(21)Lij=1+1ηij2⋅i=1,2;j=1,2.We introduce a diagonal matrix(22)L¯i=1Li100001Li200001Li100001Li2,and its inverse matrix(23)L¯i−1=Li10000Li20000Li10000Li2.(20) can be simplified into a matrix product of (15) and (22):(24)Na¯i1=1secφi⋅a¯i1⋅L¯i.The corresponding inverse matrix of (20) can be simplified using (18) and (23):(25)Na¯i1−1=1secφi⋅a¯i1⋅L¯i−1=secφi⋅L¯i−1⋅a¯i1−1.We can also write (24) and (25) in terms of the cotangent functions of rotation angles from (16) and (19):(26)Na¯i2=1cscφi⋅a¯i2⋅L¯i,Na¯i2−1=1cscφi⋅a¯i2⋅L¯i−1=cscφi⋅L¯i−1⋅a¯i2−1,


### 2.3. Reflection and Transmission Coefficients for the Anisotropic Interface in the Ice Sheet

We use the 2 × 2 reflection matrixR and the 2 × 2 transmission matrix T for the problem of EM propagation at the anisotropic interface. Considering an interface between two anisotropic ice-sheet layers, the reflection matrix R can be defined to relate the reflected upgoing waves to the downgoing incident waves. The transmission matrix T relates the transmitted downgoing waves to the downgoing incident waves. Above the interface, the total wave field consists of the incident downgoing wave and the reflected upgoing wave. Below the interface, the total wave field has only downgoing transmitted waves. The reflection and transmission at the interface of the anisotropic half-space satisfy the following relation:(27)Na¯1⋅R¯I=Na¯20T¯.We do matrix operation on equation (27) and get(28)Na¯2−1⋅Na¯1⋅R¯I=0T¯,where(29)R¯=R11R12R21R22,T¯=T11T12T21T22.Formula (28) illustrates that the eigenvector matrix of the first layer Na¯1 and the inverse eigenvector matrix of the second layer Na¯2−1 are necessary for calculating R¯ and T¯. Since the eigenvector matrix in each layer has two representations, there are four combinations to construct the matrix A¯=Na¯2−1⋅Na¯1. We first consider a form with (20) and (25):(30)A¯=Na¯21−1⋅Na¯11.Substituting (24) and (25) into (30) yields(31)A¯=Na¯21−1⋅Na¯11=secφ2⋅L¯2−1⋅a¯21−1⋅1secφ1⋅a¯11⋅L¯1=secφ2secφ1⋅L¯2−1⋅a¯21−1⋅a¯11⋅L¯1=secφ2secφ1⋅L¯2−1⋅B¯⋅L¯1.where(32)a¯21−1=12sec2φ2tgφ21η21−η21⋅tgφ21−tgφ2−η22⋅tgφ2−η22tgφ21−η21η21⋅tgφ21−tgφ2η22⋅tgφ2η22,a¯11=tgφ11tgφ111−tgφ11−tgφ11η11−1η12⋅tgφ1−1η111η12⋅tgφ1−1η11⋅tgφ1−1η121η11⋅tgφ11η12.Decomposing the eigenvector matrix into a diagonal matrix and an elementary matrix yields(33)B¯=a¯21−1⋅a¯11=Cof1⋅E¯2⋅b¯2−1⋅b¯1⋅E¯1,where(34)E¯2=120−100120−11201001201=12I−I12II,b¯2−1=2tgφ22002−2tgφ20000−η21η22⋅tgφ200η21⋅tgφ2η22,b¯1=2tgφ12002−2tgφ10000−1η111η12⋅tgφ1001η11⋅tgφ11η12,E¯1=120120012012−10100−101=12I12I−II,Cof1=12sec2φ2.I is the identity matrix, and b¯2−1 and b¯1 are block diagonal matrices. Since the product of two block diagonal matrices is still block diagonal, we arrive at(35)b¯2−1.b¯1=41+tgφ1tgφ24tgφ2−tgφ100−4tgφ2−tgφ141+tgφ1tgφ20000η21η111+tgφ1tgφ2η21η12tgφ2−tgφ100−η22η11tgφ2−tgφ1η22η121+tgφ1tgφ2.=c¯100c¯2.Therefore,(36)B¯=a¯2−1⋅a¯1=Cof1⋅E2⋅b¯2⋅b¯1⋅E1=Cof1⋅14c¯1+c¯214c¯1−c¯214c¯1−c¯214c¯1+c¯2=Cof2⋅1+η21η111+η21η12tgφ2−φ1−1+η22η11tgφ2−φ11+η22η121−η21η111−η21η12tgφ2−φ1−1−η22η11tgφ2−φ11−η22η121−η21η111−η21η12tgφ2−φ1−1−η22η11tgφ2−φ11−η22η121+η21η111+η21η12tgφ2−φ1−1+η22η11tgφ2−φ11+η22η12,=Cof2B¯aB¯bB¯bB¯a.where(37)Cof2=1+tgφ1⋅tgφ22sec2φ2.Let(38)L¯1=1L11001L12000000001L11001L12=F¯100F¯1,(39)L¯2=1L21001L22000000001L21001L22=F¯200F¯2,and(40)F¯1=1L11001L12=1/1+1η112001/1+1η122,(41)F¯2=1L21001L22=1/1+1η212001/1+1η222,One can simplify (31) by substituting (36), (38), and (39):(42)A=secφ2secφ1⋅L¯2i−1⋅B¯⋅L¯1.=secφ2secφ1.F¯2−100F¯2−1⋅Cof2B¯aB¯bB¯bB¯a⋅F¯100F¯1.=Cof3.F¯2−1⋅B¯a⋅F¯1F¯2−1⋅B¯b⋅F¯1F¯2−1⋅B¯b⋅F¯1F¯2−1⋅B¯a⋅F¯1.=Cof3⋅A¯aA¯bA¯bA¯a,where(43)Cof3=secφ2secφ1⋅Cof2=1+tgφ1tgφ22secφ1secφ2=12sinφ2−φ1.Substituting (42) into (28) yields(44)Cof3⋅A¯aA¯bA¯bA¯aRI=0T,where one can solve for the reflection and transmission coefficients R and T:(45)R¯=−A¯a−1⋅A¯b,(46)T¯=Cof3⋅A¯b⋅R¯+A¯a=Cof3⋅−A¯b⋅A¯a−1⋅A¯b+A¯a.However, the expressions ofR and T involve inverse matrix A¯a−1, which is typically difficult to obtain. Decomposing the matrix A¯ into block diagonal matrices and elementary matrices as shown in (42) and substituting them into (45) and (46) yield(47)R¯=−A¯a−1⋅A¯b=−F¯2−1⋅B¯a⋅F¯1−1⋅F¯2−1⋅B¯b⋅F¯1.=−F¯1−1⋅B¯a−1⋅F¯2.⋅F¯2−1⋅B¯b⋅F¯1=−F¯1−1⋅B¯a−1⋅B¯b⋅F¯1,(48)T¯=Cof3⋅A¯b⋅R¯+A¯a=Cof3⋅−A¯b⋅A¯a−1⋅A¯b+A¯a.=Cof3⋅F¯2−1⋅−B¯b⋅B¯a−1⋅B¯b+B¯a⋅F¯1.In this way, complex calculation ofR and T can be decomposed into product or sum between the matrix on the diagonal matrix F and the submatrix of B and its inverse matrix, where the inverse of matrix B¯a can be obtained from (36):(49)B¯a−1=1B¯a1+η22η12−1+η21η12tgφ2−φ11+η22η11tgφ2−φ11+η21η11.Leta=η21/η11,b=η21/η12,c=η22/η11,d=η22/η12, and(50)B¯a=1+η21η111+η22η12+1+η21η121+η22η11tg2φ2−φ1.=1+a1+d+1+b1+ctg2φ2−φ1,(51)R¯1=B¯a−1⋅B¯b=−1B¯a⋅1+d1−a+1+b1−ctg2φ2−φ12d−btgφ2−φ12c−atgφ2−φ11+a1−d+1+c1−btg2φ2−φ1.We can simplify the reflection coefficient (47) using (51):(52)R¯=−F¯1−1⋅R¯1⋅F¯1=−1B¯a⋅1+d1−a+1+b1−ctg2φ2−φ12d−btgφ2−φ11+1/η1121+1/η1222c−atgφ2−φ11+1/η1221+1/η1121+a1−d+1+c1−btg2φ2−φ1.Similarly, substituting(53)T¯1=−B¯b⋅B¯a−1⋅B¯b+B¯a=4B¯a⋅a1+d+a+bctg2φ2−φ1b+adtgφ2−φ1+b1+ctg3φ2−φ1−c+adtgφ2−φ1−c1+btg3φ2−φ1d1+a+d+bctg2φ2−φ1into (48) simplifies the transmission coefficient:(54)T¯=Cof3⋅F¯2−1⋅T¯1⋅F¯1=2⋅cosφ2−φ1B¯a⋅F1⋅a1+d+a+bctg2φ2−φ1F2⋅b+adtgφ2−φ1+b1+ctg3φ2−φ1F3⋅−c+adtgφ2−φ1−c1+btg3φ2−φ1F4⋅d1+a+d+bctg2φ2−φ1,where(55)F1=1+1/η2121+1/η112,F2=1+1/η2121+1/η122,F3=1+1/η2221+1/η112,F4=1+1/η2221+1/η122,B¯a=1+a1+d+1+b1+ctg2φ2−φ1.For the other three combinations ofA¯=Na¯2−1⋅Na¯1, we can obtain the similar expression for the reflection and transmission coefficients using ctgφ2−φ1. For the special case of Δφ=φ2−φ1 approximating kπ+π/2, the value of tgφ2−φ1 tends to infinity. So, we need to select appropriate combination for avoiding the infinity.(56)R¯=−1B¯a⋅1+d1−actg2φ2−φ1+1+b1−c2d−bctgφ2−φ11+1/η1121+1/η1222c−actgφ2−φ11+1/η1221+1/η1121+a1−dctg2φ2−φ1+1+c1−b,T¯=2⋅sinφ2−φ1B¯a⋅F1⋅a1+dctg3φ2−φ1+a+bc⋅ctgφ2−φ1F2⋅b+adctg2φ2−φ1+b1+cF3⋅−c+adctg2φ2−φ1−c1+bF4⋅d1+a⋅ctg3φ2−φ1+d+bcctgφ2−φ1,where(57)B¯a=1+a1+dctg2φ2−φ1+1+b1+c.


### 2.4. Special Case: Isotropic and Anisotropic Interface

We now consider the special case of isotropic ice media, where the permittivity tensor is isotropic, and (6) can be simplified to(58)Γi=000−μ000μ000ε00−ε000.Following the aforementioned derivations, the corresponding eigenvalues and eigenvectors are(59)λ=λ1λ2λ3λ4=εμ0εμ0−εμ0−εμ0,a¯1=010110101η0−1η00−1η01η,a¯2=1010010101η0−1η−1η01η0.Notice that the eigenvectors in the isotropic ice sheet only depend on the deflection angle at 0 andπ/2. Therefore, the reflection and transmission coefficients (52) and (54) can be simplified to(60)R¯=−1−a1+acos2φ2+1−d1+dsin2φ2d−a1+a1+dsin2⋅φ2d−a1+a1+dsin2⋅φ21−a1+asin2φ2+1−d1+dcos2φ2,(61)T¯=12⋅F1⋅a1+a⋅cosφ2F2⋅a1+a⋅sinφ2−F3⋅d1+d⋅sinφ2F4⋅d1+d⋅cosφ2.


## 2.1. Simplification of the Electromagnetic Model

RES survey and deep ice core analysis have proved that, in small scale, the deep ice sheet can be simplified as horizontal stratified ice layers. The ice sheet is often composed of four main fabric types [6, 27].Figure1 demonstrates the Schmidt diagrams with the four fabrics in the ice sheet. C-axis is a direction orthogonal to the crystal basal plane and often acts as the unique indicator of a crystal. The gray points or filled areas in Figure 1 indicate the possible C-axis directions or the gathering shape of C-axis (fabric) in polycrystalline ice. The symbols x and y indicate the principal axes in the horizontal plane, and the z-axis is perpendicular to the horizontal plane. Random ice (Figure 1(a)) often appears in the shallow part of the ice sheet because the stress is faint, and the C-axis points any direction kept while snow fell on the surface. The perfect single pole (Figure 1(b)) and vertical girdle (Figure 1(d)) are the two extreme states of fabric with ice crystal deformation under the stress and shear. The elongated single pole (Figure 1(c)) is the transient state between the perfect single pole and the vertical girdle under shear. Mathematically, anisotropic crystals in the ice sheet can be described by the second-order permittivity tensor, which has a diagonal form along the principal axis direction:(1)ε⟶=ε1000ε2000ε3.Figure 1
Schmidt diagrams demonstrate typical fabrics in the ice sheet: (a) a random ice, (b) a perfect single pole, (c) an elongated single pole, and (d) a vertical girdle.
(a)(b)(c)(d)The fabric can be categorized based on the signatures of differentiation of permittivity along principal axes. The permittivity of a crystal along or perpendicular to itsC-axis is denoted as ε⊥ and ε//, respectively, which had precisely been measured in the laboratory [20, 21]. The four main fabric types can be categorized through the permittivity tensor: a random ice (Figure 1(a)) has permittivity elements of ε1=ε2=ε3; a perfect single-pole ice (Figure 1(b)) has permittivity elements of ε1=ε2=ε⊥ and ε3=ε//; a vertical girdle fabric (Figure 1(d)) can be achieved when the permittivity elements satisfy ε1=ε// and ε2=ε3=ε//+ε⊥/2; and the permittivity tensor of an elongated single pole (Figure 1(c)) has ε1<ε2≪ε3.In radar analysis, EM wave propagations are typically described in two coordinates: measurement and media coordinates. Measurement coordinate is used to describe the amplitude and propagation direction of EM components. Media coordinate is used to define the principal axes of permittivity tensor and geometry of media. Figure2 presents the relative position of the measurement and media coordinates, where x′ denotes the semimajor axis of the indicatrix ellipses and y′ denotes the semiminor axis. Note that the z-axis in the measurement coordinate coincides with the z’-axis in the media coordinate because EM waves propagate as plane waves with normal incident angle into the ice sheet, where isochrones are generally horizontal. In the layered ice model, the measurement coordinate is shared among all layers, whereas the media coordinate differs in each layer and is simply a rotation of the measurement coordinate with an azimuth angle φaroundz/z’axis.Figure 2
Schematic diagram demonstrates the measurement coordinate ({x, y, z}) and media coordinate ({x′,y′,z′}). They own a common z/z' axis and rotate with azimuth angle φ. The red ellipse is the indicatrix ellipse, which indicates the principal axis direction of the permittivity tensor. Here, we suppose that the semimajor axis denotes x′.Applying the coordinate transformation to the diagonal permittivity tensor (formula (1)) yields the generic form of the permittivity tensor in the measurement coordinate:(2)ε¯=ε1sin2φ+ε2cos2φ12ε1−ε2sin2φ012ε1−ε2sin2φε1cos2φ+ε2sin2φ000ε3.


## 2.2. Eigenanalysis of EM Wave Propagations in Anisotropic Ice

The EM wave propagation in anisotropic ice-sheet layers can be described by differential-form Maxwell’s equations in time-harmonic field:(3)∇×E⇀i=−jωμ0H⇀i,∇×H⇀i=jωε¯iE⇀i,where E⇀i and H⇀i are the electric and magnetic fields in the i-th ice sheet layer, ω is the angular frequency, μ0 is the isotropic permeability of the i-th layer, and ε¯i is the permittivity tensor of the i-th layer. In the measurement coordinate, we suppose z-axis is the direction of EM propagation. It satisfies the case of normal incidence. So, electric field E and magnetic field H just exist in the plane of x- and y-axis determining. Combining electric and magnetic fields into one vector(4)ψi=ExEyHxHyTyields the matrix form of Maxwell’s equations [43–45]:(5)dψizdz=jω⋅Γiz⋅ψiz,where(6)Γi=000−μ000μ00ε1−ε2sinφcosφε1cos2φ+ε2sin2φ00−ε1sin2φ+ε2cos2φ−ε1−ε2sinφcosφ00.The entries ofΓi are expressed in terms of the azimuth angle φ and permittivity tensor in the i-th layer. For vertically propagating plane waves with wave number λ, the eigensolutions of (5) can be found by letting(7)ψ=ψ0eλz,which yields an eigenvalue problem with the eigenvalues λ:(8)Γi−λI⋅ψ0=0.In general, linear system (8) has four eigenvalues. Therefore, the general solution of (7) can be written as the linear combination of eigenvectors of (8):(9)ψz=A1a1ejλ1z+A2a2ejλ2z+A3a3e−jλ3z+A4a4e−jλ4z,or equivalently in matrix-vector notation:(10)ψz=a¯ejλ¯zA¯,where(11)a¯=a1a2a3a4,A¯=A1A2A3A4T,ejλ¯z=ejλ1zejλ2zejλ3zejλ4z.We can further reformulate (8) into(12)Γ−λI⋅X=0,where(13)λ⟶=λ1λ2λ3λ4=ε1μ0ε2μ0−ε1μ0−ε2μ0.Equation (12) has two fundamental solutions:(14)XI=x1x2x3x4=1ctgφλiμ0ctgφ−λiμ0,XII=x1x2x3x4=tanφ1λiμ0−λiμ0tanφ.For the normal incident EM wave along the -z/z′ axis, the wave field can be decomposed into waves I and II, where waves I and II are in the direction of x′ and y′ axis, respectively. In the measurement coordinate, waves I and wave II have rotation angles of φ and φ+π/2, respectively. The eigenvectors of the wave field in the i-th layer have two forms:(15)a¯i1=tgφi1tgφi11ctgπ2+φi1ctgπ2+φiεi1μ0εi2μ0⋅ctgπ2+φi−εi1μ0−εi2μ0⋅ctgπ2+φi−εi1μ0⋅tgφi−εi2μ0εi1μ0⋅tgφiεi2μ0,=tgφi1tgφi11−tgφi1−tgφi1ηi1−1ηi2⋅tgφi−1ηi11ηi2⋅tgφi−1ηi1⋅tgφi−1ηi21ηi1⋅tgφi1ηi2.(16)a¯i2=1tgπ2+φi1tgπ2+φictgφi1ctgφi1εi1μ0⋅ctgφiεi2μ0−εi1μ0⋅ctgφi−εi2μ0−εi1μ0−εi2μ0⋅tgπ2+φiεi1μ0εi2μ0⋅tgπ2+φi.=1−ctgφi1−ctgφictgφi1ctgφi11ηi1⋅ctgφi1ηi2−1ηi1⋅ctgφi−1ηi2−1ηi11ηi2⋅ctgφi1ηi1−1ηi2⋅ctgφi,where(17)ηi1=μ0εi1,ηi2=μ0εi2,are intrinsic impedances of waves I and II in the i-th layer, subscripts 1 and 2 represent waves I and II, respectively, and φi is the azimuthal rotation angle of the i-th layer in the measurement coordinate with respect to the media coordinate. Subscripts (1) and (2) indicate the two forms of matrices in (15) and (16), which are eigenvectors composed of tangent and cotangent functions of the azimuth angles, respectively. Taking the inverse of matrices (15) and (16) yields(18)a¯i1−1=12sec2φitgφi1tgφi11−tgφi1−tgφiηi1−ηi2⋅tgφi−ηi1ηi2⋅tgφi−ηi1⋅tgφi−ηi2ηi1⋅tgφiηi2T,(19)a¯i2−1=12csc2φi1−ctgφi1−ctgφictgφi1ctgφi1ηi1⋅ctgφiηi2−ηi1⋅ctgφi−ηi2−ηi1ηi2⋅ctgφiηi1−ηi2⋅ctgφiT.The inverse matrices (18) and (19) resemble the matrices (15) and (16), except for a transpose operator, inverse of wave impedance, and a scaling factor. Therefore, the inverse matrices (18) and (19) can be easily written once the dielectric tensor and azimuth angle are given.The calculation of reflection-transmission coefficients needs normalized eigenvectors, so we first normalize each column of the matrix in (15) by the L2 norm of the column vector:(20)Na¯i1=1secφitgφLi11Li2tgφiLi11Li21Li1−tgφiLi21Li1−tgφiLi21ηi1⋅Li1−tgφiηi2⋅Li2−1ηi1⋅Li1tgφiηi2⋅Li2−tgφηi1⋅Li1−1ηi2⋅Li2tgφiηi1⋅Li11ηi2⋅Li2,where(21)Lij=1+1ηij2⋅i=1,2;j=1,2.We introduce a diagonal matrix(22)L¯i=1Li100001Li200001Li100001Li2,and its inverse matrix(23)L¯i−1=Li10000Li20000Li10000Li2.(20) can be simplified into a matrix product of (15) and (22):(24)Na¯i1=1secφi⋅a¯i1⋅L¯i.The corresponding inverse matrix of (20) can be simplified using (18) and (23):(25)Na¯i1−1=1secφi⋅a¯i1⋅L¯i−1=secφi⋅L¯i−1⋅a¯i1−1.We can also write (24) and (25) in terms of the cotangent functions of rotation angles from (16) and (19):(26)Na¯i2=1cscφi⋅a¯i2⋅L¯i,Na¯i2−1=1cscφi⋅a¯i2⋅L¯i−1=cscφi⋅L¯i−1⋅a¯i2−1,


## 2.3. Reflection and Transmission Coefficients for the Anisotropic Interface in the Ice Sheet

We use the 2 × 2 reflection matrixR and the 2 × 2 transmission matrix T for the problem of EM propagation at the anisotropic interface. Considering an interface between two anisotropic ice-sheet layers, the reflection matrix R can be defined to relate the reflected upgoing waves to the downgoing incident waves. The transmission matrix T relates the transmitted downgoing waves to the downgoing incident waves. Above the interface, the total wave field consists of the incident downgoing wave and the reflected upgoing wave. Below the interface, the total wave field has only downgoing transmitted waves. The reflection and transmission at the interface of the anisotropic half-space satisfy the following relation:(27)Na¯1⋅R¯I=Na¯20T¯.We do matrix operation on equation (27) and get(28)Na¯2−1⋅Na¯1⋅R¯I=0T¯,where(29)R¯=R11R12R21R22,T¯=T11T12T21T22.Formula (28) illustrates that the eigenvector matrix of the first layer Na¯1 and the inverse eigenvector matrix of the second layer Na¯2−1 are necessary for calculating R¯ and T¯. Since the eigenvector matrix in each layer has two representations, there are four combinations to construct the matrix A¯=Na¯2−1⋅Na¯1. We first consider a form with (20) and (25):(30)A¯=Na¯21−1⋅Na¯11.Substituting (24) and (25) into (30) yields(31)A¯=Na¯21−1⋅Na¯11=secφ2⋅L¯2−1⋅a¯21−1⋅1secφ1⋅a¯11⋅L¯1=secφ2secφ1⋅L¯2−1⋅a¯21−1⋅a¯11⋅L¯1=secφ2secφ1⋅L¯2−1⋅B¯⋅L¯1.where(32)a¯21−1=12sec2φ2tgφ21η21−η21⋅tgφ21−tgφ2−η22⋅tgφ2−η22tgφ21−η21η21⋅tgφ21−tgφ2η22⋅tgφ2η22,a¯11=tgφ11tgφ111−tgφ11−tgφ11η11−1η12⋅tgφ1−1η111η12⋅tgφ1−1η11⋅tgφ1−1η121η11⋅tgφ11η12.Decomposing the eigenvector matrix into a diagonal matrix and an elementary matrix yields(33)B¯=a¯21−1⋅a¯11=Cof1⋅E¯2⋅b¯2−1⋅b¯1⋅E¯1,where(34)E¯2=120−100120−11201001201=12I−I12II,b¯2−1=2tgφ22002−2tgφ20000−η21η22⋅tgφ200η21⋅tgφ2η22,b¯1=2tgφ12002−2tgφ10000−1η111η12⋅tgφ1001η11⋅tgφ11η12,E¯1=120120012012−10100−101=12I12I−II,Cof1=12sec2φ2.I is the identity matrix, and b¯2−1 and b¯1 are block diagonal matrices. Since the product of two block diagonal matrices is still block diagonal, we arrive at(35)b¯2−1.b¯1=41+tgφ1tgφ24tgφ2−tgφ100−4tgφ2−tgφ141+tgφ1tgφ20000η21η111+tgφ1tgφ2η21η12tgφ2−tgφ100−η22η11tgφ2−tgφ1η22η121+tgφ1tgφ2.=c¯100c¯2.Therefore,(36)B¯=a¯2−1⋅a¯1=Cof1⋅E2⋅b¯2⋅b¯1⋅E1=Cof1⋅14c¯1+c¯214c¯1−c¯214c¯1−c¯214c¯1+c¯2=Cof2⋅1+η21η111+η21η12tgφ2−φ1−1+η22η11tgφ2−φ11+η22η121−η21η111−η21η12tgφ2−φ1−1−η22η11tgφ2−φ11−η22η121−η21η111−η21η12tgφ2−φ1−1−η22η11tgφ2−φ11−η22η121+η21η111+η21η12tgφ2−φ1−1+η22η11tgφ2−φ11+η22η12,=Cof2B¯aB¯bB¯bB¯a.where(37)Cof2=1+tgφ1⋅tgφ22sec2φ2.Let(38)L¯1=1L11001L12000000001L11001L12=F¯100F¯1,(39)L¯2=1L21001L22000000001L21001L22=F¯200F¯2,and(40)F¯1=1L11001L12=1/1+1η112001/1+1η122,(41)F¯2=1L21001L22=1/1+1η212001/1+1η222,One can simplify (31) by substituting (36), (38), and (39):(42)A=secφ2secφ1⋅L¯2i−1⋅B¯⋅L¯1.=secφ2secφ1.F¯2−100F¯2−1⋅Cof2B¯aB¯bB¯bB¯a⋅F¯100F¯1.=Cof3.F¯2−1⋅B¯a⋅F¯1F¯2−1⋅B¯b⋅F¯1F¯2−1⋅B¯b⋅F¯1F¯2−1⋅B¯a⋅F¯1.=Cof3⋅A¯aA¯bA¯bA¯a,where(43)Cof3=secφ2secφ1⋅Cof2=1+tgφ1tgφ22secφ1secφ2=12sinφ2−φ1.Substituting (42) into (28) yields(44)Cof3⋅A¯aA¯bA¯bA¯aRI=0T,where one can solve for the reflection and transmission coefficients R and T:(45)R¯=−A¯a−1⋅A¯b,(46)T¯=Cof3⋅A¯b⋅R¯+A¯a=Cof3⋅−A¯b⋅A¯a−1⋅A¯b+A¯a.However, the expressions ofR and T involve inverse matrix A¯a−1, which is typically difficult to obtain. Decomposing the matrix A¯ into block diagonal matrices and elementary matrices as shown in (42) and substituting them into (45) and (46) yield(47)R¯=−A¯a−1⋅A¯b=−F¯2−1⋅B¯a⋅F¯1−1⋅F¯2−1⋅B¯b⋅F¯1.=−F¯1−1⋅B¯a−1⋅F¯2.⋅F¯2−1⋅B¯b⋅F¯1=−F¯1−1⋅B¯a−1⋅B¯b⋅F¯1,(48)T¯=Cof3⋅A¯b⋅R¯+A¯a=Cof3⋅−A¯b⋅A¯a−1⋅A¯b+A¯a.=Cof3⋅F¯2−1⋅−B¯b⋅B¯a−1⋅B¯b+B¯a⋅F¯1.In this way, complex calculation ofR and T can be decomposed into product or sum between the matrix on the diagonal matrix F and the submatrix of B and its inverse matrix, where the inverse of matrix B¯a can be obtained from (36):(49)B¯a−1=1B¯a1+η22η12−1+η21η12tgφ2−φ11+η22η11tgφ2−φ11+η21η11.Leta=η21/η11,b=η21/η12,c=η22/η11,d=η22/η12, and(50)B¯a=1+η21η111+η22η12+1+η21η121+η22η11tg2φ2−φ1.=1+a1+d+1+b1+ctg2φ2−φ1,(51)R¯1=B¯a−1⋅B¯b=−1B¯a⋅1+d1−a+1+b1−ctg2φ2−φ12d−btgφ2−φ12c−atgφ2−φ11+a1−d+1+c1−btg2φ2−φ1.We can simplify the reflection coefficient (47) using (51):(52)R¯=−F¯1−1⋅R¯1⋅F¯1=−1B¯a⋅1+d1−a+1+b1−ctg2φ2−φ12d−btgφ2−φ11+1/η1121+1/η1222c−atgφ2−φ11+1/η1221+1/η1121+a1−d+1+c1−btg2φ2−φ1.Similarly, substituting(53)T¯1=−B¯b⋅B¯a−1⋅B¯b+B¯a=4B¯a⋅a1+d+a+bctg2φ2−φ1b+adtgφ2−φ1+b1+ctg3φ2−φ1−c+adtgφ2−φ1−c1+btg3φ2−φ1d1+a+d+bctg2φ2−φ1into (48) simplifies the transmission coefficient:(54)T¯=Cof3⋅F¯2−1⋅T¯1⋅F¯1=2⋅cosφ2−φ1B¯a⋅F1⋅a1+d+a+bctg2φ2−φ1F2⋅b+adtgφ2−φ1+b1+ctg3φ2−φ1F3⋅−c+adtgφ2−φ1−c1+btg3φ2−φ1F4⋅d1+a+d+bctg2φ2−φ1,where(55)F1=1+1/η2121+1/η112,F2=1+1/η2121+1/η122,F3=1+1/η2221+1/η112,F4=1+1/η2221+1/η122,B¯a=1+a1+d+1+b1+ctg2φ2−φ1.For the other three combinations ofA¯=Na¯2−1⋅Na¯1, we can obtain the similar expression for the reflection and transmission coefficients using ctgφ2−φ1. For the special case of Δφ=φ2−φ1 approximating kπ+π/2, the value of tgφ2−φ1 tends to infinity. So, we need to select appropriate combination for avoiding the infinity.(56)R¯=−1B¯a⋅1+d1−actg2φ2−φ1+1+b1−c2d−bctgφ2−φ11+1/η1121+1/η1222c−actgφ2−φ11+1/η1221+1/η1121+a1−dctg2φ2−φ1+1+c1−b,T¯=2⋅sinφ2−φ1B¯a⋅F1⋅a1+dctg3φ2−φ1+a+bc⋅ctgφ2−φ1F2⋅b+adctg2φ2−φ1+b1+cF3⋅−c+adctg2φ2−φ1−c1+bF4⋅d1+a⋅ctg3φ2−φ1+d+bcctgφ2−φ1,where(57)B¯a=1+a1+dctg2φ2−φ1+1+b1+c.


## 2.4. Special Case: Isotropic and Anisotropic Interface

We now consider the special case of isotropic ice media, where the permittivity tensor is isotropic, and (6) can be simplified to(58)Γi=000−μ000μ000ε00−ε000.Following the aforementioned derivations, the corresponding eigenvalues and eigenvectors are(59)λ=λ1λ2λ3λ4=εμ0εμ0−εμ0−εμ0,a¯1=010110101η0−1η00−1η01η,a¯2=1010010101η0−1η−1η01η0.Notice that the eigenvectors in the isotropic ice sheet only depend on the deflection angle at 0 andπ/2. Therefore, the reflection and transmission coefficients (52) and (54) can be simplified to(60)R¯=−1−a1+acos2φ2+1−d1+dsin2φ2d−a1+a1+dsin2⋅φ2d−a1+a1+dsin2⋅φ21−a1+asin2φ2+1−d1+dcos2φ2,(61)T¯=12⋅F1⋅a1+a⋅cosφ2F2⋅a1+a⋅sinφ2−F3⋅d1+d⋅sinφ2F4⋅d1+d⋅cosφ2.


## 3. Numerical Results

It had been mentioned that echo difference exists in radar exploration while rotating the antenna horizontally [19, 22–24]. Therefore, we consider two numerical models to verify our aforementioned formulae for computing the reflection and transmission coefficients. The aforementioned formulae have no special limit to the range of permittivity. However, for test and verification of the real reflections in the ice sheet, we use the permittivities ε⊥=3.152 and ε//=3.189 for model parameters, which had precisely been measured in the laboratory [20, 21].


### 3.1. Reflection and Transmission Coefficients for an Interface between Anisotropic Layers

Firstly, we construct a two-layer ice model, where both layers are anisotropic and the values of the entry permittivity tensors areε1=3.152, ε2=3.189, and ε3=3.189. In model 1, we fix media coordinates for both layers and keep Δφ=φ2−φ1=π/6, where φ1 and φ2 are azimuth angles of the first and second layer relative to the measurement coordinate. Then, we rotate the measurement coordinate from 0 to 2π and calculate each reflection and transmission coefficient at each rotating angle. This operation emulates the polarizing antenna rotates around the horizontal plane inversely.The reflection and transmission coefficients are drawn in Figure3. Both reflection and transmission coefficients are kept constant, while rotating the measurement coordinate from 0 to 2π in the horizontal plane. When we change the permittivity or geometry Δφ, the values of the reflection and transmission coefficients will change too, but they still remain constant in a new level while measurement coordinate rotating.Figure 3
The reflection and transmission coefficients from an interface between two homogeneous anisotropic ice layers. The measurement coordinate rotates from 0 to2π in the horizontal plane.
(a)(b)In model 2, we fix the measurement coordinate and scan the difference of media azimuth angleΔφ=φ2−φ1 from 0 to 2π, where φ1 and φ2 are azimuth angles of the first and second layer. In this case, we observe that the reflection and transmission coefficients vary periodically in Figure 4.Figure 4
The reflection and transmission coefficients in the anisotropy-anisotropy interface. The upper layer remains fixed, and the lower layer rotates from 0 to 2π in the horizontal plane.
(a)(b)


### 3.2. Reflection and Transmission Coefficients for an Interface between Isotropic and Anisotropic Layers

Considering the real survey environment of the ice sheet, we construct another two-layer ice model, in which the upper layer is isotropic withε1=ε2=ε3=1 as the free air layer, and the lower layer is anisotropic with ε1=3.152, ε2=3.189, and ε3=3.189 as the real ice layer. In this model, we fix media coordinates for both layers and keep φ2−φ1=0 and φ2−φ1=−π/6 (Figure 5), where φ1 and φ2 are azimuth angles of ice fabric in the first and second layer. We rotate the measurement coordinate from 0 to 2π in the horizontal plane.Figure 5
The reflection and transmission coefficients at the interface between free air and anisotropic ice layers are calculated. The difference of azimuth angles in the upper and lower layer is set atφ2−φ1=0 and φ2−φ1=−π/6, respectively. The measurement coordinate rotates from 0 to 2π.
(a)(b)(c)(d)In this case, we observe that reflection and transmission coefficients vary periodically in Figure5 while we rotate the measurement coordinates. Comparing the variation of reflection and transmission, we can see the curves shift right −π/6, while the difference of azimuth angle changes from 0 to −π/6. This can be used for tracking the azimuth of the maximum echo in radar survey on the ice sheet.


## 3.1. Reflection and Transmission Coefficients for an Interface between Anisotropic Layers

Firstly, we construct a two-layer ice model, where both layers are anisotropic and the values of the entry permittivity tensors areε1=3.152, ε2=3.189, and ε3=3.189. In model 1, we fix media coordinates for both layers and keep Δφ=φ2−φ1=π/6, where φ1 and φ2 are azimuth angles of the first and second layer relative to the measurement coordinate. Then, we rotate the measurement coordinate from 0 to 2π and calculate each reflection and transmission coefficient at each rotating angle. This operation emulates the polarizing antenna rotates around the horizontal plane inversely.The reflection and transmission coefficients are drawn in Figure3. Both reflection and transmission coefficients are kept constant, while rotating the measurement coordinate from 0 to 2π in the horizontal plane. When we change the permittivity or geometry Δφ, the values of the reflection and transmission coefficients will change too, but they still remain constant in a new level while measurement coordinate rotating.Figure 3
The reflection and transmission coefficients from an interface between two homogeneous anisotropic ice layers. The measurement coordinate rotates from 0 to2π in the horizontal plane.
(a)(b)In model 2, we fix the measurement coordinate and scan the difference of media azimuth angleΔφ=φ2−φ1 from 0 to 2π, where φ1 and φ2 are azimuth angles of the first and second layer. In this case, we observe that the reflection and transmission coefficients vary periodically in Figure 4.Figure 4
The reflection and transmission coefficients in the anisotropy-anisotropy interface. The upper layer remains fixed, and the lower layer rotates from 0 to 2π in the horizontal plane.
(a)(b)


## 3.2. Reflection and Transmission Coefficients for an Interface between Isotropic and Anisotropic Layers

Considering the real survey environment of the ice sheet, we construct another two-layer ice model, in which the upper layer is isotropic withε1=ε2=ε3=1 as the free air layer, and the lower layer is anisotropic with ε1=3.152, ε2=3.189, and ε3=3.189 as the real ice layer. In this model, we fix media coordinates for both layers and keep φ2−φ1=0 and φ2−φ1=−π/6 (Figure 5), where φ1 and φ2 are azimuth angles of ice fabric in the first and second layer. We rotate the measurement coordinate from 0 to 2π in the horizontal plane.Figure 5
The reflection and transmission coefficients at the interface between free air and anisotropic ice layers are calculated. The difference of azimuth angles in the upper and lower layer is set atφ2−φ1=0 and φ2−φ1=−π/6, respectively. The measurement coordinate rotates from 0 to 2π.
(a)(b)(c)(d)In this case, we observe that reflection and transmission coefficients vary periodically in Figure5 while we rotate the measurement coordinates. Comparing the variation of reflection and transmission, we can see the curves shift right −π/6, while the difference of azimuth angle changes from 0 to −π/6. This can be used for tracking the azimuth of the maximum echo in radar survey on the ice sheet.


## 4. Conclusion and Discussion

### 4.1. Physical Meaning of Eigenvalues and Eigenvectors

The eigenvalues and eigenvectors contain information about the relation between electromagnetic components and physical properties of media when the EM wave propagates in anisotropic media. Multiplying (13) by the angular frequency ω yields the upgoing and downgoing wave vectors of I and II waves:(62)ω⋅λ=ω⋅λ1λ2λ3λ4=ωε1μ0ωε2μ0−ωε1μ0−ωε2μ0=kIkII−kI−kII.The column vectors of eigenvector matrices in (15) and (16) show the proportional relation of EM components, ExEyHxHy. The four columns in eigenvector matrices represent the normalized component of the upgoing and downgoing wave vectors corresponding to I and II waves, respectively. In anisotropic media,(63)ExiEyi=tgφ,HxiHyi=tgφ,i=I,II(64)ExiHyi=ηi,EyiHxi=ηi.i=I,IIEquation (63) shows that the ratios Ex/Ey and Hx/Hy are only functions of the media angle φ, which is measured between the measurement and the media coordinate. The intrinsic impedance of media relates the magnetic and electric fields, as shown in (64). So, the eigenvector matrix describes the proportion relationship between EM components of I and II waves and proves that eigenvectors are only determined by intrinsic impedance of media and the deflection angle φ. Therefore, we can conclude that the eigenvector matrix describes the proportional relation of EM components of I and II waves. Moreover, the eigenvectors are only determined by the intrinsic impedance of media and the deflection angle φ.


### 4.2. Simplification of the Formulae for Reflection and Transmission Coefficients

For natural media,F¯1 and F¯2 in formulas (40) and (41) tend to the unit matrix. In formula (52),(65)1+1/η1121+1/η122≈1.In formula (54), F1, F2, F3, and F4 also tend to 1, so formulas (52) and (54) can be simplified to(66)R¯=−F¯1−1⋅R¯1⋅F¯1=−1B¯a⋅1+d1−a+1+b1−ctg2φ2−φ12d−btgφ2−φ12c−atgφ2−φ11+a1−d+1+c1−btg2φ2−φ1,T¯=Cof3⋅F¯2−1⋅T¯1⋅F¯1=2⋅cosφ2−φ1B¯a⋅a1+d+a+bctg2φ2−φ1b+adtgφ2−φ1+b1+ctg3φ2−φ1−c+adtgφ2−φ1−c1+btg3φ2−φ1d1+a+d+bctg2φ2−φ1,B¯a=1+a1+d+1+b1+ctg2φ2−φ1.


### 4.3. Polarimetric Radar Survey and Geometry Analysis

A typical configuration of polarimetric radar consists of two orthogonal transmitters and two orthogonal receivers. Four transmitter-receiver pairs are possible, either copolarized or cross-polarized, which are shown in Figure6. Conventional polarimetric radar implements configuration with time-share single transmitter to double receiver alternately. Orthogonal echo differences can be used to infer the anisotropic feature of the media.Figure 6
The schematic of polarimetric radar about the four transmitter-receiver pairs.For an interface between isotropic-anisotropic layers, we can assume that the measurement coordinate is {X, Y} in Figure 7, and the media coordinate of the anisotropic layer has deflection angle φ2 relative to the measurement coordinate. EH is the incident wave in the isotropic layer. It can be decomposed into EHI and EHII along the semimajor and semiminor axes of the indicatrix ellipse. These two incident waves generate two reflection waves ERHI and ERHII from the interface:(67)RI=E⟶HIRE⟶HI=η2I−η1η2I+η1=ε1−ε2Iε1+ε2I,RII=E⟶HIIRE⟶HII=η2II−η1η2II+η1=ε1−ε2IIε1+ε2II.Figure 7
The electric field components of the incident wave and reflection wave in the isotropic and anisotropic interface. Blue ellipse represents the indicatrix ellipse, and its long axis points to the direction of I; wave and short axis points to ; wave. Just incident wave ofEH is chosen to demonstrate here, incident wave EV in axis y can be decomposed with the same method.Therefore, the reflection from the interface along thex-axis is the vector sum of ERHI and ERHII along the x-axis. The reflection coefficient R11 can be calculated as(68)R11=E⟶HRE⟶H=E⟶HIxR+E⟶HIIxRE⟶H=E⟶HIR⋅cosφ2+E⟶HIIR⋅sinφ2E⟶H=RI⋅E⟶H⋅cos2φ2+RII⋅E⟶H⋅sin2φ2E⟶H=RI⋅cos2φ2+RII⋅sin2φ2,and similarly, the reflection coefficient R12 can be computed as(69)R12=E⟶HIyR+E⟶HIIyRE⟶H=E⟶HIR⋅sinφ2−E⟶HIIR⋅cosφ2E⟶H=RI⋅E⟶H⋅sinφ2⋅cosφ2+RII⋅E⟶H⋅sinφ2⋅cosφ2E⟶H=12⋅sin2⋅φ2RI−RII.For incident waveEV along the y-axis, we can derive the reflection coefficients R21 and R22 in a similar fashion, which have been shown in (60).


## 4.1. Physical Meaning of Eigenvalues and Eigenvectors

The eigenvalues and eigenvectors contain information about the relation between electromagnetic components and physical properties of media when the EM wave propagates in anisotropic media. Multiplying (13) by the angular frequency ω yields the upgoing and downgoing wave vectors of I and II waves:(62)ω⋅λ=ω⋅λ1λ2λ3λ4=ωε1μ0ωε2μ0−ωε1μ0−ωε2μ0=kIkII−kI−kII.The column vectors of eigenvector matrices in (15) and (16) show the proportional relation of EM components, ExEyHxHy. The four columns in eigenvector matrices represent the normalized component of the upgoing and downgoing wave vectors corresponding to I and II waves, respectively. In anisotropic media,(63)ExiEyi=tgφ,HxiHyi=tgφ,i=I,II(64)ExiHyi=ηi,EyiHxi=ηi.i=I,IIEquation (63) shows that the ratios Ex/Ey and Hx/Hy are only functions of the media angle φ, which is measured between the measurement and the media coordinate. The intrinsic impedance of media relates the magnetic and electric fields, as shown in (64). So, the eigenvector matrix describes the proportion relationship between EM components of I and II waves and proves that eigenvectors are only determined by intrinsic impedance of media and the deflection angle φ. Therefore, we can conclude that the eigenvector matrix describes the proportional relation of EM components of I and II waves. Moreover, the eigenvectors are only determined by the intrinsic impedance of media and the deflection angle φ.


## 4.2. Simplification of the Formulae for Reflection and Transmission Coefficients

For natural media,F¯1 and F¯2 in formulas (40) and (41) tend to the unit matrix. In formula (52),(65)1+1/η1121+1/η122≈1.In formula (54), F1, F2, F3, and F4 also tend to 1, so formulas (52) and (54) can be simplified to(66)R¯=−F¯1−1⋅R¯1⋅F¯1=−1B¯a⋅1+d1−a+1+b1−ctg2φ2−φ12d−btgφ2−φ12c−atgφ2−φ11+a1−d+1+c1−btg2φ2−φ1,T¯=Cof3⋅F¯2−1⋅T¯1⋅F¯1=2⋅cosφ2−φ1B¯a⋅a1+d+a+bctg2φ2−φ1b+adtgφ2−φ1+b1+ctg3φ2−φ1−c+adtgφ2−φ1−c1+btg3φ2−φ1d1+a+d+bctg2φ2−φ1,B¯a=1+a1+d+1+b1+ctg2φ2−φ1.


## 4.3. Polarimetric Radar Survey and Geometry Analysis

A typical configuration of polarimetric radar consists of two orthogonal transmitters and two orthogonal receivers. Four transmitter-receiver pairs are possible, either copolarized or cross-polarized, which are shown in Figure6. Conventional polarimetric radar implements configuration with time-share single transmitter to double receiver alternately. Orthogonal echo differences can be used to infer the anisotropic feature of the media.Figure 6
The schematic of polarimetric radar about the four transmitter-receiver pairs.For an interface between isotropic-anisotropic layers, we can assume that the measurement coordinate is {X, Y} in Figure 7, and the media coordinate of the anisotropic layer has deflection angle φ2 relative to the measurement coordinate. EH is the incident wave in the isotropic layer. It can be decomposed into EHI and EHII along the semimajor and semiminor axes of the indicatrix ellipse. These two incident waves generate two reflection waves ERHI and ERHII from the interface:(67)RI=E⟶HIRE⟶HI=η2I−η1η2I+η1=ε1−ε2Iε1+ε2I,RII=E⟶HIIRE⟶HII=η2II−η1η2II+η1=ε1−ε2IIε1+ε2II.Figure 7
The electric field components of the incident wave and reflection wave in the isotropic and anisotropic interface. Blue ellipse represents the indicatrix ellipse, and its long axis points to the direction of I; wave and short axis points to ; wave. Just incident wave ofEH is chosen to demonstrate here, incident wave EV in axis y can be decomposed with the same method.Therefore, the reflection from the interface along thex-axis is the vector sum of ERHI and ERHII along the x-axis. The reflection coefficient R11 can be calculated as(68)R11=E⟶HRE⟶H=E⟶HIxR+E⟶HIIxRE⟶H=E⟶HIR⋅cosφ2+E⟶HIIR⋅sinφ2E⟶H=RI⋅E⟶H⋅cos2φ2+RII⋅E⟶H⋅sin2φ2E⟶H=RI⋅cos2φ2+RII⋅sin2φ2,and similarly, the reflection coefficient R12 can be computed as(69)R12=E⟶HIyR+E⟶HIIyRE⟶H=E⟶HIR⋅sinφ2−E⟶HIIR⋅cosφ2E⟶H=RI⋅E⟶H⋅sinφ2⋅cosφ2+RII⋅E⟶H⋅sinφ2⋅cosφ2E⟶H=12⋅sin2⋅φ2RI−RII.For incident waveEV along the y-axis, we can derive the reflection coefficients R21 and R22 in a similar fashion, which have been shown in (60).



---

*Source: 9869873-2020-06-05.xml*