Multi-channel and multi-polarization radar measurements around the NEEM site

Abstract

Ice properties inferred from multi-polarization measurements, such as birefringence and crystal orientation fabric (COF), can provide insight into ice strain, viscosity, and ice flow. In 2008, the Center for Remote Sensing of Ice Sheets (CReSIS) used a ground-based VHF (very high frequency) radar to take multi-channel and multi-polarization measurements around the NEEM (North Greenland Eemian Ice Drilling) site. The system operated with 30 MHz bandwidth at a center frequency of 150 MHz. This paper describes the radar system, antenna configurations, data collection, and processing and analysis of this data set. Within the framework derived from uniaxial ice crystal model, we found that ice birefringence dominates the power variation patterns of co-polarization and cross-polarization measurements in the area of 100 km${}^{\mathrm{2}}$ around the ice core site. The phase shift between ordinary and extraordinary waves increases nonlinearly with depth. The ice optic axis lies in planes that are close to the vertical plane and perpendicular or parallel to the ice divide depending on depth. The ice optic axis has an average tilt angle of about $\mathrm{11.6}{}^{\circ }$ vertically, and its plane may rotate either clockwise or counterclockwise by about $\mathrm{10}{}^{\circ }$ across the 100 km${}^{\mathrm{2}}$ area, and at a specific location the plane may rotate slightly counterclockwise as depth increases. Comparisons between the radar observations, simulations, and ice core fabric data are in very good agreement. We calculated the effective colatitude at different depths by using azimuth and colatitude measurements of the $c$ axis of ice crystals. We obtained an average effective $c$ axis tilt angle of $\mathrm{9.6}{}^{\circ }$ from the vertical axis, very comparable to the average optic axis tilt angle estimated from radar polarization measurements. The comparisons give us confidence in applying this polarimetric radio echo sounding technique to infer profiles of ice fabric in locations where there are no ice core measurements.

Introduction

Radio echo sounding of ice sheets is a decades-old technique (Robin et al., 1969). Ice thickness profiles, internal layers, bed topography, and basal conditions revealed by radar echo sounding are widely used for analysis and modeling of ice mass balance, temperature, flow, and dynamics (Shepherd et al., 2012; MacGregor, 2015; Bamber et al., 2013; Hughes et al., 2016). Many previous studies show that ice properties, such as the birefringence and crystal orientation fabric (COF), could be inferred from multi-polarization radar measurements. The strength and polarization of radar echoes are affected by ice COF with respect to antenna orientations (Campbell and Orange, 1974; Bentley, 1975). Ice COF is closely related to ice flow history and has a strong effect on ice deformation rate (Azuma and Higashi, 1985); it is one of three reasons for radio echo reflections from internal ice layers (Fujita and Mae, 1994; Fujita et al., 1999; Eisen et al., 2007). Hargreaves (1977) first proposed the theory of polarization and propagation of electromagnetic waves in a birefringent medium. He developed a method for determining radar wave polarizations by rotating a pair of orthogonal receiving antennas and observing variations in echo strength, and presented a linear incremental trend in the phase difference between ordinary and extraordinary waves caused by the ice birefringence from field measurements in the vicinity of DYE-3 in Greenland. Thereafter, many researchers further investigated this subject (Woodruff and Doake, 1979; Doake, 1981; Doake et al., 2002; Fujita and Mae, 1994; Fujita et al., 1999, 2006; Matsuoka et al., 2003, 2009, 2012; Drews et al., 2012). Woodruff and Doake (1979) suggested that large changes in the birefringence of ice sheets at the grounding zone were caused by the effect of tidal strain on crystal orientation. This technique could distinguish the floating from the grounded areas of the ice sheets. Fujita and Mae (1994) proposed a dual-frequency technique to distinguish permittivity-based and conductivity-based reflections from internal ice layers. Matsuoka et al. (2003) made multi-polarization radar measurements at two frequencies (60 and 179 MHz) from Dome Fuji toward the coast in East Antarctica, inferring the relationship between the identified COF spatial patterns of anisotropic reflectivity and birefringence to ice flow dynamics. Drews et al. (2012) tested and discussed in detail birefringence and other types of anisotropy in natural ice (e.g., anisotropic distribution of elongated air inclusions) for their effect on radar data by applying different scattering models.

The NEEM project (North Greenland Eemian Ice Drilling), led by Denmark, was an international collaboration of 14 nations between 2007 and 2011 with the purpose of retrieving an ice core from the previous Eemian interglacial to understand the dynamics of past climates under conditions similar to the warming climate we are approaching. In the past, these ice core data have been used to infer the ice sheet's response to the strong warming in the early Eemian (NEEM community members, 2013) and derive ice and gas chronology (Rasmussen et al., 2013). The fabric measured from the NEEM ice core and comparisons with GRIP (Greenland Ice Core Project) and NGRIP (North Greenland Ice Core Project) ice cores may help to constrain ice flow modeling along the ridge from GRIP to NEEM and onwards to Camp Century (Weikusat and Kipfstuhl, 2010; Eichler et al., 2013; Montagnat et al., 2014). In 2007, CReSIS (the Center for Remote Sensing of Ice Sheets) conducted both airborne and ground-based radar surveys to help identify the best area for the NEEM site. In 2008, CReSIS used an improved ground-based VHF radar capable of mapping the deepest layers in surveys along gridded lines covering about 100 km${}^{\mathrm{2}}$ around the NEEM site (Blake et al., 2009). Multi-channel and polarimetric measurements were collected with special antenna configurations along circular paths to infer changes in the ice fabric from radio echo sounding. The polarization measurements of this data set were preliminarily analyzed in 2015 (Vélez González, 2015). In this paper, we present further analysis to show that the ice birefringence dominates power variations in the received echoes from internal layers as the antennas rotate along the circular paths. We find that the phase shift between the ordinary and extraordinary waves increases nonlinearly with depth because of the ice birefringence. We infer the ice optic axis at the NEEM site is close to the vertical planes, approximately perpendicular or parallel to the ice divide for different depth ranges, and may rotate within a range between $-\mathrm{10}{}^{\circ }$ (clockwise) and $+\mathrm{10}{}^{\circ }$ (counterclockwise) across the survey area. At a specific location, the plane also rotates slightly counterclockwise as the depth increases.

This paper first describes the surface-based radar system and antenna configurations used in the data collection of multi-polarization measurements, and gives details of the data processing techniques and analysis. Next, it compares the measurements with simulation results according to the theory of plane wave propagation in a birefringent medium. Lastly, this study provides a comparison between radio echo sounding measurements and ice core observations with discussions and conclusions on the ice COF and birefringence properties around the NEEM ice core site.

Multi-polarization measurements

Radar system overview and antenna configurations

The radar system used in the survey was a modified ground-based version of the multi-channel radar depth sounder (MCRDS) developed by CReSIS in 2006 and 2007 (Lohoefener, 2006; Marathe, 2008). Table 1 summarizes the key parameters of the radar system. There is a detailed description of the system found in Supplement Sect. S1 alongside a simplified system block diagram in Fig. S1.

System parameters.

Parameter	Value	Units
Frequency band	135 to 165	MHz
Pulse duration	3 or 10	$\mathrm{\mu }$s
Pulse repetition frequency	10	kHz
Peak transmit power	400	W
Receiver noise figure	3.9	dB
Loop sensitivity	212	dB
Minimum detectable signal	$-\mathrm{162}$	dBm
Range resolution	2.8	m
Coherent averages	64
Antenna element	log-periodic
Antenna gain	6.7	dB
A/D converter	12	bit
Sampling rate	120	MSPS

Two antenna configurations were used during field operations: HH polarization setup for non-polarimetric measurements and quad-polarimetric setup for multi-polarization measurements (Fig. 1). In the HH polarization setup, both transmit and receive antennas were H polarized (i.e., the antenna plane is in the along-track direction). The left and right transmitters (TX1 and TX2) were in front of the eight evenly spaced receivers (RX1–RX8). Each transmitter consisted of a pair of log-periodic antennas fed simultaneously by the RF transmitter to form a common phase center. In the quad-polarimetric setup, the antenna of TX2 and RX5–RX8 were V polarized by rotating them 90${}^{\circ }$ about the vertical axis (Fig. 1, right). We installed the antennas on a custom-made sled towed by a tracked vehicle (Fig. S3, left). We also installed a Topcon Legacy E GPS receiver between the transmitters and receivers to record geolocations with NMEA $GPGGA messages every 0.1 s. We provide details about antenna geometry, radiation patterns, and field installation in Supplement Sect. S2.

Data collection

We retrieved multi-channel, multi-waveform, and multi-polarization data with the quad-polarimetric antenna setup on 16 August 2008. Figure 2a illustrates data collection paths. The survey paths began by the NEEM ice core camp and went southwest perpendicular to the ice divide, passed through turning points 1, 2, 3, and 4 in sequence and returned to the starting point. The section from point 2 to 3 is along the ice divide at an angle of $-\mathrm{60.7}$${}^{\circ }$ from north. The distances between locations 1 and 4 and between locations 2 and 3 are about 10 km. At each turning point, we traveled in a circular path with a radius of about 30 m to take quad-polarimetric measurements at all orientations (Fig. 2b). The circular path at location 2 was not closed as at the other three locations.

During the data collection, the two transmitters TX1 and TX2 alternatively transmitted four waveforms (3 $\mathrm{\mu }$s pulse by TX2, 3 $\mathrm{\mu }$s pulse by TX1, 10 $\mathrm{\mu }$s pulse by TX2 and 10 $\mathrm{\mu }$s pulse by TX1) in sequence at a PRF (pulse repetition frequency) of 10 kHz. The signals from the eight receive channels were digitized separately. The 64 consecutive pulses were stacked together and saved as one record. We thus took multi-channel measurements of four polarizations (HH, HV, VH, and VV) simultaneously. The speed of the tracked vehicle was about 3 m s${}^{-\mathrm{1}}$, resulting in a distance of about 0.077 m between two data records.

Antenna configurations: (a) non-polarimetric HH polarization setup; (b) quad-polarimetric setup.

[Figure omitted. See PDF]

Paths of multi-polarization measurements: (a) paths along ice divide and perpendicular to ice divide. (b) Circular paths at turns. The paths are plotted in the local ENU (east–north–up) coordinate system at the NEEM ice core camp site (77.45${}^{\circ }$ N, 51.06${}^{\circ }$ W), marked by a red dot at ($\mathrm{0},\mathrm{0}$). The travel direction is illustrated by the arrows along the paths. The local ice divide is along the path from location 2 to 3, as marked by the dashed line with arrows on both ends. The red dot and the dashed line with arrows in the inset map show the location of the NEEM site and the direction of the local ice divide in Greenland. The cross formed by a red and a blue segment shows the location of a crossover where the data were collected with the non-polarimetric HH polarization antenna configuration (see discussions on Fig. 7).

[Figure omitted. See PDF]

Data processing and analysis

Receive channel equalization

While most polarization measurements previously reported in the literature were single channel data, we collected multi-channel polarization measurements for our purposes in this research. Although we can choose data from any one of the channels for analysis, we perform our analysis using combined multi-channel data because it has higher SNR (signal-to-noise ratio), reduced clutter and speckle artifacts, and clearly revealed deep ice layers. For a radar system with multiple receive channels and antenna arrays, any difference in hardware (e.g., cable length, amplifiers, attenuators, and antenna matching/feed network) may result in channel mismatches in two-way propagation time delay, amplitude, and phase. Channel equalization – the process to compensate these channel mismatches – is required to combine the multichannel data to increase the signal-to-noise ratio and reduce the across-track clutter. Supplement Sect. S3 presents the receive channel equalization method and results.

Compared to single-channel measurements, the combined multi-channel measurements with channel equalization have more distinct internal layers, less speckle, and higher SNR. As mentioned earlier, these attributes are especially critical for the analysis of deep layers.

Transmit power mismatch

As discussed in Supplement Sect. S4, we identified a 12.7 dB transmit power mismatch between TX1 and TX2. Although the transmit power towards nadir by TX2 was reduced, the system had adequate performance for polarimetric measurements because of the redundancy of the quad-polarimetric configuration. Hereinafter, we will only analyze polarization measurements from transmit TX1, which includes HH and HV transmit–receive polarizations.

Multi-polarization measurements along Circle 3: (a) HH, transmit TX1, receive channels RX1, RX2, RX3, and RX4 combined; (b) HV, transmit TX1, receive channels RX5, RX6, RX7, and RX8 combined. The horizontal axis is the angle in degrees with respect to the local north, the vertical axis is ice depth in kilometers, and the color bar represents the relative power level in decibels.

[Figure omitted. See PDF]

Circular data analysis

Figure 3 presents the HH and HV polarization measurements at Circle 3. For HH measurements, the data from RX1, RX2, RX3, and RX4 were averaged together after compensating for inner-channel mismatches (as discussed in Sect. 3.1). The same procedure was used for HV measurements with data from RX5, RX6, RX7, and RX8. To further increase the signal-to-noise ratio and to reduce speckle, two data records were coherently averaged along the circle, followed by four incoherent integrations. The coherent and incoherent integrations resulted in a sampling interval of about 0.54 m along the circle, corresponding to an azimuth change of about $\mathrm{1.125}{}^{\circ }$.

In Fig. 3, the top $\sim \mathrm{127}$ m is blank. This corresponds to the eclipsing time, which is equal to the duration of the 3 $\mathrm{\mu }$s chirps that the radar has to wait before recording data for unambiguous range detection. The sudden intensity drop between the depth of 1450 m and the ice bed at 2500 m corresponds to the so-called echo free zone (EFZ). As seen in Fig. 3, the receive power level obtained for HH measurements is higher than that of HV measurements for shallow depths between 127 and 750 m, and the power level of HV measurements becomes comparable to that of HH measurements for depths of 750 to 1450 m. Even in the EFZ, we observe three distinct layers at depths of 1753, 1832, and 1894 m, respectively. The most unique and impressive feature in Fig. 3 is the obvious signal extinction patterns separated by about $\mathrm{90}{}^{\circ }$, in both the HH and HV measurements of middle depth range (between 750 and 1450 m).

Next, we will compare results from polarimetric measurements with simulations based on existing theory using the uniaxial ice crystal model (Hargreaves, 1977) and show that the above patterns are caused by dominant ice birefringence effects. We will treat the nadir-transmitted radar signals into the ice sheets as plane waves traveling downward. If we take the nadir as the positive $z$ direction, then the electric field of the plane wave $\mathit{E}$ is in the $x$–$y$ plane perpendicular to the $z$ axis. Its components $E_x$ and $E_y$ are given as $$\begin{array}{}& {\displaystyle }& {\displaystyle }{E}_x=A_x\cos (\mathit{\omega }t-kz),& {\displaystyle }& {\displaystyle }{E}_y=A_y\cos (\mathit{\omega }t-kz+\mathrm{\varnothing }),\end{array}$$ where $A_x$ and $A_y$ are the amplitude of the $x$ component and $y$ component, respectively, $\mathit{\omega }$ is the angular frequency, $k$ is the wave number and $\mathrm{\varnothing }$ is the phase difference between the two components. $\mathit{E}$ is a function of time at any fixed position $z$, and the vector tip will trace a curve in the $x$–$y$ plane as time elapses. The plane wave is of linear, circular, or elliptical polarization when this curve is, respectively, a straight line, a circle, or an ellipse (Ulaby, 1981). In general, Eqs. (1) and (2) describe an elliptical polarization. Circular polarization is a special case when$A_x=A_y\ne \mathrm{0}$ and $\mathrm{\varnothing }=\pm \mathrm{90}{}^{\circ }$ and linear polarization is a special case when $\mathrm{\varnothing }=\mathrm{0}$ and $A_x$ and $A_y$ are not zeros, or when either $A_x$ or $A_y$ is zero.

The usual dipole and Yagi antennas radiate linear polarized waves (in general, $\mathrm{\varnothing }=\mathrm{0}$ and $A_x$, and $A_y$ are not zeros). The polarization and antenna planes match. Two identical and orthogonally crossed dipoles, fed with a $\mathrm{90}{}^{\circ }$ phase difference, radiate circular polarized waves (Balanis, 2005). The 17-element log-periodic antennas radiate linear polarized waves and the polarization planes match with the antenna planes. Ideally, if there is no polarization change of radio waves during propagation, there is no signal loss from polarization mismatch if the transmit (TX) and receive (RX) antennas are co-polarized (HH or VV); otherwise, there is a signal loss because of the polarization mismatch between the TX and RX antennas, with cross-polarized TX–RX as the extreme case (HV or VH). Previous studies have found that the polarization state of radio waves might change after propagating through and reflected from ice sheets because of ice properties such as the birefringence of ice (Bogorodsky et al., 1976). In these cases, the signal loss would occur from polarization mismatch, even if the TX and RX antennas are co-polarized.

Ice sheets are polycrystalline formed by many ice crystals of different shapes and sizes at different orientations. The orientation of an ice crystal can be determined by its $c$ axis orientation. The crystal orientation fabrics (COFs) in ice sheets are classified based on the $c$ axis projection distribution on Schmidt diagrams. Vertical single-maximum, elongated single-maximum, and vertical girdle are three typical COF types usually found in ice cores from ice sheets, mainly drilled at divides and domes (Alley et al., 1997; Matsuoka et al., 2003; Fujita et al., 2006; Eisen et al., 2007; Montagnat et al., 2012; Weikusat et al., 2017; see Faria et al., 2014a for a review). For single-maximum COF, all $c$ axes align up in a cone centered along or close to the $z$ axis, and in the classical glaciological projection into the horizontal the points in Schmidt diagrams (stereographic projections) concentrate close to the center of the circle as this center represents the vertical direction ($z$ axis). The point distribution of elongated single-maximum COF is an ellipse around the center of the circular diagram. A band across the diagram is called vertical girdle COF, representing a great-circle in the unprojected reality. These three COF cases correspond to thinning dominated flow or radial dome-flow situations (vertical single maximum COF), diverging or parallel flow (girdle COF), and general flow situations with diverging and thinning components (elongated maximum COF) in ice sheets, as ice crystals rotate differently under these various compression and extension situations (Faria et al., 2014b; Paterson, 1994).

Different COF types result in different dielectric permittivity properties which can be described by the dielectric permittivity tensor ${\in }_{\mathrm{p}}$ in the principle coordinate system: $${\in }_{\mathrm{p}}=\left[\begin{array}{ccc}{\mathit{ϵ}}_{x_{\mathrm{p}}}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& {\mathit{ϵ}}_{y_{\mathrm{p}}}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& {\mathit{ϵ}}_{z_{\mathrm{p}}}\end{array}\right],$$ where $z_{\mathrm{p}}$ is parallel with the symmetric axis of the COF, ${\mathit{ϵ}}_{x_{\mathrm{p}}}$ and ${\mathit{ϵ}}_{y_{\mathrm{p}}}$ are the two components in the $x_{\mathrm{p}}$–$y_{\mathrm{p}}$ plane that is perpendicular to $z_{\mathrm{p}}$. For single-maximum COF, ${\mathit{ϵ}}_{x_{\mathrm{p}}}\approx {\mathit{ϵ}}_{y_{\mathrm{p}}}$; for elongated single-maximum COF, $\left|{\mathit{ϵ}}_{z_{\mathrm{p}}}-{\mathit{ϵ}}_{x_{\mathrm{p}}}\right|$ or $\left|{\mathit{ϵ}}_{z_{\mathrm{p}}}-{\mathit{ϵ}}_{y_{\mathrm{p}}}\right|\gg \left|{\mathit{ϵ}}_{x_{\mathrm{p}}}-{\mathit{ϵ}}_{y_{\mathrm{p}}}\right|>\mathrm{0}$; and for vertical-girdle COF, $\left|{\mathit{ϵ}}_{z_{\mathrm{p}}}-{\mathit{ϵ}}_{x_{\mathrm{p}}}\right|$ or $\left|{\mathit{ϵ}}_{z_{\mathrm{p}}}-{\mathit{ϵ}}_{y_{\mathrm{p}}}\right|>\left|{\mathit{ϵ}}_{x_{\mathrm{p}}}-{\mathit{ϵ}}_{y_{\mathrm{p}}}\right|>\mathrm{0}$. For ice with elongated single-maximum COF, the dielectric permittivity is about the same around the $z_{\mathrm{p}}$ axis and the difference between ${\mathit{ϵ}}_{z_{\mathrm{p}}}$ and ${\mathit{ϵ}}_{x_{\mathrm{p}}}$ or ${\mathit{ϵ}}_{y_{\mathrm{p}}}$ is about 0.034 according to previous study (Matsuoka et al., 1997). The anisotropic property in the ice dielectric permittivity tensor is known as the birefringence of ice sheets.

For convenience, if the ice is birefringent, we can choose the $x$ axis and $y$ axis so that the optic axis of the ice is in the $y$–$z$ plane and along the principal axis $z_{\mathrm{p}}$, which is at an angle of $\mathit{\beta }$ to the $z$ axis (see Fig. 12c). Under the assumption that the ice is uniaxial about $z_{\mathrm{p}}$ (i.e., ${\mathit{ϵ}}_{x_{\mathrm{p}}}={\mathit{ϵ}}_{y_{\mathrm{p}}})$, according to the solutions of Maxwell's equations, a plane wave transmitted vertically down the $z$ axis to a birefringent ice sheet would split into an ordinary wave (its electric field is in the direction of the $x$ axis) and an extraordinary wave (its electric field is in the $y$–$z$ plane). The wave numbers of the two waves $k_{\mathrm{1}}$ and $k_{\mathrm{2}}$ would be different, which results in a phase difference at depth $z$ given by $$\mathrm{\varnothing }={\mathrm{\varnothing }}_{\mathrm{2}}-{\mathrm{\varnothing }}_{\mathrm{1}}=\mathrm{2}(k_{\mathrm{2}}-k_{\mathrm{1}})z,$$ where $k_{\mathrm{1}}$ and $k_{\mathrm{2}}$ are given by Eq. (A4a) and (A4b) in Hargreaves, (1977). Therefore, for $\mathit{\beta }\ne \mathrm{0}$, and the transmit polarization not aligned with the $x$ or $y$ axis, the transmitted wave of linear polarization would become elliptically polarized after passing through the birefringent ice because of the above phase shift.

When denoting the angle between the TX1 plane and the $x$ axis as $\mathit{\alpha }$, the received reflections from ice layers for a receiver positioned at an angle of $\mathit{\gamma }$ with respect to the $x$ axis would be $$\begin{array}{}& {\displaystyle }& {\displaystyle }{E}_x^{\prime }=A_x^{\prime }\cos \mathit{\alpha }\cos \mathit{\gamma }\cos (\mathit{\omega }t-kz),& {\displaystyle }& {\displaystyle }{E}_y^{\prime }=A_y^{\prime }\sin \mathit{\alpha }\sin \mathit{\gamma }\cos (\mathit{\omega }t-kz+\mathrm{\varnothing }),\end{array}$$ where $A_x^{\prime }$ and $A_y^{\prime }$ are the amplitude of the $x$ and $y$ components, respectively, which include the losses from propagation, ice attenuation, and reflections. Therefore, the power of the received signals is $$\begin{array}{cc}{\displaystyle }{P}_z& {\displaystyle }\left(\mathit{\alpha },\mathit{\theta },\mathrm{\varnothing }\right)=E{{}^{\prime }}_x^{\mathrm{2}}+E{{}^{\prime }}_y^{\mathrm{2}}=A{{}^{\prime }}_x^{\mathrm{2}}{\cos }^{\mathrm{2}}\mathit{\alpha }{\cos }^{\mathrm{2}}\mathit{\gamma }\\ & {\displaystyle }& {\displaystyle }+A{{}^{\prime }}_y^{\mathrm{2}}{\sin }^{\mathrm{2}}\mathit{\alpha }{\sin }^{\mathrm{2}}\mathit{\gamma }+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}A{{}^{\prime }}_{x}A{{}^{\prime }}_y\sin \mathrm{2}\mathit{\alpha }\sin \mathrm{2}\mathit{\gamma }\cos \mathrm{\varnothing },\end{array}$$ For RX1–RX4, the antenna plane is parallel with that of TX1 (co-polarization), so $\mathit{\gamma }=\mathit{\alpha }$; for RX5–RX8, the antenna plane is perpendicular to that of TX1 (cross-polarization), so $\mathit{\gamma }=\mathit{\alpha }+\mathrm{90}{}^{\circ }$. Figure 4a presents the simulation results of the received signal power as a function of $\mathit{\alpha }$ and $\mathrm{\varnothing }$ for the co-polarization and cross-polarization cases, respectively, based on Eq. (7) assuming $A_x^{\prime }=A_y^{\prime }=\mathrm{1}$. The patterns repeat every $\mathrm{90}{}^{\circ }$ of $\mathit{\alpha }$ and are thus displayed only in the range of 0–90${}^{\circ }$ for $\mathit{\alpha }$. Each pattern is symmetric with respect to $\mathrm{45}$ and $\mathrm{180}{}^{\circ }$ for $\mathit{\alpha }$ and $\mathrm{\varnothing }$, respectively. The co-polarization power has minima at $\mathit{\alpha }=\mathrm{45}{}^{\circ }$ and $\mathrm{\varnothing }=\mathrm{180}{}^{\circ }$, and the cross-polarization power's maxima are at the same locations.

Effects of the angle $\mathit{\alpha }$ between TX1 and $x$ axis and birefringence phase shift $\mathrm{\varnothing }$ on isotropic and anisotropic reflected power. The color bars in panels (a) and (b) give the power scale in the range from 0 to $-\mathrm{25}$ dB. In panel (c), the blue and red lines correspond to co-polarization and cross-polarization, respectively; the solid lines are simulated isotropic pattern, dashed lines are simulated anisotropic pattern, and dotted lines are filtered measurements along Circle 3.

[Figure omitted. See PDF]

The solid lines in Fig. 4c compare the simulated power level of the received signals of the co-polarization and cross-polarization for phase shifts of 30, 60, 90, and 120${}^{\circ }$ (the phase shift increases with depth) as $\mathit{\alpha }$ rotates from 0 to 360${}^{\circ }$. The extinctions of cross-polarization power occur at $\mathit{\alpha }=\mathrm{0}$, 90, 180, 270, and 360${}^{\circ }$. At $\mathrm{\varnothing }=\mathrm{30}$${}^{\circ }$, the co-polarization power level varies little as $\mathit{\alpha }$ rotates. The minima at $\mathit{\alpha }=n\mathrm{90}{}^{\circ }+\mathrm{45}{}^{\circ }(n=\mathrm{0},\mathrm{1},\mathrm{2},\mathrm{3})$ are only 0.3 dB below the maxima (0 dB) at $\mathit{\alpha }=n\mathrm{90}{}^{\circ }(n=\mathrm{0},\mathrm{1},\mathrm{2},\mathrm{3},\mathrm{4})$; the maxima of the cross-polarization power at $\mathit{\alpha }=n\mathrm{90}{}^{\circ }+\mathrm{45}{}^{\circ }(n=\mathrm{0},\mathrm{1},\mathrm{2},\mathrm{3})$ are 11.75 dB down with respect to the maximum co-polarization power. At $\mathrm{\varnothing }=\mathrm{60}{}^{\circ }$, the co-polarization power variations become visible as $\mathit{\alpha }$ rotates. The minima are 1.25 dB down from the maxima, and the maxima of the cross-polarization power increase to $-\mathrm{6}$ dB. At $\mathrm{\varnothing }=\mathrm{90}{}^{\circ }$, the co-polarization power has minima of $-\mathrm{3}$ dB, and the maxima of cross-polarization power match the co-polarization minima. At $\mathrm{\varnothing }=\mathrm{120}{}^{\circ }$, the co-polarization power has minima of $-\mathrm{6}$ dB and the maxima of the cross-polarization power increase to $-\mathrm{1.25}$ dB and become comparable with the maximum co-polarization power.

Power profiles with respect to the angle of north at different depths for HH and HV measurements. Four colors (blue, red, green, and magenta) are cycled through for distinguishing close profiles and all the profiles are numbered from 1 to 14 in increasing order, according to their respective depths from the top to the bottom. The annotations in panel (a) illustrate the isotropic and anisotropic patterns discussed with Fig. 7 at the end of Sect. 3.3.

[Figure omitted. See PDF]

Previous surveys involving polarimetric measurements were performed along a short straight line by moving back and forth with fixed antenna orientation, and repeating the measurements after changing the antenna orientation (e.g., Fujita et al., 1999; Matsuoka et al., 2003). In this way, multiple measurements along a straight line for a fixed antenna orientation could be stacked together to reduce the fading of received echoes. In our case, we cannot directly stack measurements together along the circular paths because the antenna orientations are continuously changing during the measurements. Therefore, we filtered the HH and HV measurements of Circle 3 with a windowed 2-D moving average filter, as described in Supplement Sect. S5.

Instead of considering the reflected power from a few specific and isolated layers at different depths (as in previous studies), we looked at the bulk properties of the received power at all depths. For every 10 range bins between 100 and 650 (corresponding to depth increments of about 28 m between depths of 276 and 1826 m), we calculated the co-polarization and cross-polarization power profiles against the angle between TX1 and the north. Figure 5 presents the HH and HV power profiles at depth increments of about 112 m from the top to the bottom of the plots. These profiles were obtained after filtering the data by Hanning windows with $M=\mathrm{51}$ and $N=\mathrm{11}$ in Supplement Eqs. (S6) and (S7).

By integrating Eq. (7) over $[\mathrm{0},\mathrm{2}\mathit{\pi }]$ of $\mathit{\alpha }$ for co-polarization and cross-polarization cases, the birefringence phase shift as a function of ice depth $z$ can be solved as $$\begin{array}{cc}{\displaystyle }\mathrm{\varnothing }=& {\displaystyle }{\cos }^{-\mathrm{1}}\left({\displaystyle \frac{\mathrm{1}-\mathrm{3}{P}_R}{\mathrm{1}+P_R}}\right),\\ & {\displaystyle }=& {\displaystyle }{\cos }^{-\mathrm{1}}\left\{\left[\mathrm{1}-{\displaystyle \frac{\frac{\mathrm{3}}{\mathrm{2}}max\left(P_{z_{\mathrm{cr}}}\right)}{{\overline{P}}_{z_{\mathrm{co}}}}}\right]/\phantom{\left[\mathrm{1}-{\displaystyle \frac{\frac{\mathrm{3}}{\mathrm{2}}max\left(P_{z_{\mathrm{cr}}}\right)}{{\overline{P}}_{z_{\mathrm{co}}}}}\right]\left[\mathrm{1}+{\displaystyle \frac{\frac{\mathrm{1}}{\mathrm{2}}max\left(P_{z_{\mathrm{cr}}}\right)}{{\overline{P}}_{z_{\mathrm{co}}}}}\right]}\left[\mathrm{1}+{\displaystyle \frac{\frac{\mathrm{1}}{\mathrm{2}}max\left(P_{z_{\mathrm{cr}}}\right)}{{\overline{P}}_{z_{\mathrm{co}}}}}\right]\right\},\end{array}$$ where $P_R={\overline{P}}_{z_{\mathrm{cr}}}/{\overline{P}}_{z_{\mathrm{co}}}$ is the ratio of the received average cross-polarization power (${\overline{P}}_{z_{\mathrm{cr}}}$) to average co-polarization power (${\overline{P}}_{z_{\mathrm{co}}}$) at depth $z$ (Eq. 9b in Hargreaves, 1977). Based on Eq. (8), within a range of depth, we could infer the birefringence phase shift using the obtained HH and HV power profiles at different depths. We ignored the power profiles at depths shallower than 324 m because the filtering uses truncated data at depths less than 254 m for the 3 $\mathrm{\mu }$s waveform. We also ignored the power profiles deeper than 925 m, in which the maximum cross-polarization power becomes larger than twice the average co-polarization power, resulting in no real solutions from Eq. (8). The red curve in Fig. 6a presents the estimated birefringence phase shift at Circle 3 as a function of depth. The plot shows that the phase shift between the ordinary wave and extraordinary wave increases nonlinearly. Between the depths of 333 and 755 m, the phase shift increases relatively slowly by $\mathrm{36}{}^{\circ }$ (from $\mathrm{50}$ to $\mathrm{86}{}^{\circ })$, and between the depths of 755 and 925 m, the phase shift increases by $\mathrm{71}{}^{\circ }$ (from $\mathrm{86}$ to $\mathrm{157}{}^{\circ }$) – almost 5 times faster. The red circles in the plot correspond to phase shifts of $\mathrm{60}$, $\mathrm{90}$, and $\mathrm{120}{}^{\circ }$ at depths of 405, 700, and 806 m, respectively. The HH and HV power profiles (blue and red dots in Fig. 4c) at these three depths are overlapped over the simulated results (blue and red solid lines). For comparison, we converted the angle between TX1 and the north in Fig. 5 to $\mathit{\alpha }$ – the angle between Tx1 and the ice divide ($\mathit{\alpha }=\mathrm{0}{}^{\circ }$ when the plane of TX1 is along the ice divide) by circularly shifting the curves to the right by $\mathrm{60.7}{}^{\circ }$ (the angle between the ice divide and the north). We also shifted the level of the HH and HV power profiles up by 93.6, 103.3, and 111.7 dB – the values of the maximum HH power below 0 dB at each depth. We observe that the variation patterns of both the HH and HV power profiles match well with the simulated ones. The HV power level is lower than the HH at $\mathrm{\varnothing }=\mathrm{60}{}^{\circ }$, becoming comparable with HH power levels at $\mathrm{\varnothing }=\mathrm{90}$ and $\mathrm{120}{}^{\circ }$. The angles where the maxima and minima occur align with the simulated ones, although there are minor offsets. Therefore, we conclude that the HH and HV power variations are mainly determined by the ice birefringence. We notice that the major mismatch between the measured power profiles and the simulated ones are with the HH, where its maxima at $\mathrm{\varnothing }=\mathrm{0},\mathrm{180}$ and $\mathrm{360}{}^{\circ }$ are relatively smaller than the ones at $\mathrm{\varnothing }=\mathrm{90}$ and $\mathrm{270}{}^{\circ }$ by $\sim \mathrm{4.5}$ dB. This mismatch comes from the anisotropic nature of the reflections, not taken into account in the simulation. To confirm the anisotropic reflection effects, we let $A_x^{\prime }=A_y^{\prime }/\mathrm{1.7}=\mathrm{1}/\mathrm{1.7}$ in Eq. (7) (corresponding to $\sim \mathrm{4.5}$ dB less power reflection along $x$ axis), and the simulated patterns are plotted in Fig. 4b to compare with the isotropic patterns in Fig. 4a. The simulated anisotropic HH and HV power profiles for phase shifts of $\mathrm{30}$, $\mathrm{60}$, $\mathrm{90}$, and $\mathrm{120}{}^{\circ }$ are plotted in Fig. 4c with blue and red dashed lines, respectively. The measured anisotropic patterns in HH power profiles match well with the simulated ones. Later in this section, we will further discuss the depth and spatial variations of the anisotropic patterns observed in HH measurements.

(a) Estimated birefringence phase shift as a function of depth; (b) HV power minimum locations. The ice divide is parallel to $\mathrm{119.3}{}^{\circ }$ and $\mathrm{299.3}{}^{\circ }$ in this figure.

[Figure omitted. See PDF]

(a) Power variation patterns and their depth profiles for the HH polarimetric measurements at the four circles. ISO indicates no obvious anisotropy observed. ANISO-1–ANISO-2 are anisotropic patterns and indicate the ice optic axis is perpendicular or parallel to the ice divide. (b) Stacked normalized power–depth profiles and their differences at a crossover from non-polarimetric measurements with HH antenna configuration. The normalization is done by dividing both profiles by the maximum of the red one at 327 m. The plot of the differences is area-filled. (c) Paths of the crossover (see the cross mark in Fig. 2a for its location relative to the four circles). (d) Anisotropic pattern ANISO-1 at the crossover. (e) Anisotropic pattern ANISO-2 at the crossover.

[Figure omitted. See PDF]

We performed similar analyses of the data obtained from measurements performed at the other three circles. Circle 2 was not completed during data collection; therefore, we used the measurements from the half circle and completed the circle by repeating the measurement according to the power variation periodicity. The derived birefringence phase shift-depth profiles are also shown in Fig. 6a in blue, green, and black curves for Circle 1, 2, and 4, respectively. The minima of the cross-polarization measurement power profiles are the most distinct and unambiguous feature shown in Figs. 3b and 5b; therefore, we tracked the antenna orientation angles between TX1 and the north, where the minima occur in four circle locations at depths not greater than 1404 m (see Fig. 6b). We observe the following major characteristics from depth and spatial variations in Fig. 6.

• The nonlinearities of birefringence phase shifts are similar for all four circles (e.g., plateaus between depths of 700 and 800 m, sudden and faster increases in the phase shift after the plateau).

• While the phase shift profiles of Circle 2, 3, and 4 are close to each other, Circle 1 is departed from the others; the phase shift of the plateau is at least 26${}^{\circ }$ higher. In Fig. 6a, we observe a faster phase shift slope before and after the plateau at Circle 1, which may suggest greater anisotropy at Circle 1 than other circle locations.

• The locations of the HV power minima at Circle 1 and 2 are very close to each other. Their offsets from the ones at Circle 3 are less than $\mathrm{10}{}^{\circ }$ (rotate clockwise) on average for most depths, and the minima at Circle 4 are constantly smaller than those at Circle 3 by around $\mathrm{10}{}^{\circ }$ (rotated counterclockwise) on average for all depths (see Fig. 12b).

• At all circles, we observe the general trend of the angle to the north of the HV minima rotating slightly counterclockwise as the depth increases, which implies that the plane of the ice optic axis rotates slightly counterclockwise as the depth increases.

  The three curves at the bottom of Fig. 5 correspond to depths of 1516, 1629, and 1742 m in the EFZ. The HV power variations still follow the birefringence patterns, and the rotation of the plane of the ice optic axis continues.

  Based on Fig. 6b, we could only infer that the ice optic axis is in a vertical plane that is either parallel or perpendicular to the ice divide. However, we could remove this ambiguity by investigating the anisotropy in the power variation patterns of HH measurements. Figure 5a marks three distinct power variation patterns of the HH measurements we identified and used to characterize the ice anisotropic properties: (a) for ISO, the ice is isotropic around axis $z$ with no obvious power peaks at shallow ice depth when birefringence phase shift is small or with four power peaks at about the same level; (b) for ANISO-1, the ice is anisotropic between $x$ and $y$ axes, with two higher power peaks and two lower peaks interleaved; (c) ANISO-2 is similar to ANISO-1 with an orientation shift of $\mathrm{90}{}^{\circ }$. Both ANISO-1 and ANISO-2 result from the effects of anisotropic reflections at internal ice layers. Because the reflection coefficient along the direction of the extraordinary wave field ($y$ axis) is larger than along the ordinary wave field ($x$ axis) and the two higher power peaks are at orientations perpendicular to the ice divide in ANISO-1, ANISO-1 indicates that the ice optic axis is in the $y$–$z$ plane perpendicular to the ice divide and ANISO-2 indicates that the ice optic axis is in the $y$–$z$ plane parallel to the ice divide.

  Figure 7a presents the depth profiles of the above three power variation patterns traced from the HH polarimetric measurements at the locations of the four circles. The transition patterns between the boundaries of any two patterns are not considered and plotted in the charts. In Fig. 7a, we can observe the following characteristics.

• 1.

  For shallow depth above $\sim \mathrm{450}$ m, no anisotropic property is observed between the $x$ and $y$ axes, and the ice is in ISO pattern.

• 2.

  Between $\sim \mathrm{450}$ and 1065 m, at all circle locations the ice is anisotropic between the $x$ and $y$ axes, and the ice optic axis is in a plane that is approximately perpendicular to the ice divide. This anisotropic property continues to $\sim \mathrm{1319}$ m at Circle 2. Below the first ANISO-1, the ice returns to ISO pattern for depths from 1094 to 1178 m, from 1263 to 1319 m, from 1094 to 1319 m, and from 1065 to 1375 m at Circle 1 to Circle 4, respectively. There is a second ANISO-1 pattern at Circle 1 and Circle 2, between 1206 and 1347 m and between 1347 and 1375 m, respectively.

• 3.

  There is an ANISO-2 pattern at all circle locations and the optic axis is parallel to the ice divide; this occurs at depths from 1375 to 1545 m, from 1404 to 1545 m, from 1347 to 1407 m and from 1404 to 1545 m at Circle 1 to Circle 4, respectively; the pattern switch from ANISO-1 to ANISO-2 at Circle 1 and Circle 2, from ISO to ANISO-2 at Circle 3 and 4.

• 4.

  For depths below ANISO-2, no anisotropic property is observed between the $x$ and $y$ axes, and the ice is in ISO pattern.

We also observed the above isotropic and anisotropic patterns in the data collected with the non-polarimetric antenna configuration (see Fig. 1a) at crossovers along grid lines that are either parallel or perpendicular to the ice divide. Figure 7b presents an example of stacked and normalized power–depth profiles at one of these crossovers, illustrating the isotropic and anisotropic patterns at different depth ranges. At this crossover, from top to bottom the ice column has ISO (324–450 m), ANISO-1 (450–1000 m), ISO (1000–1340 m), ANISO-2 (1340–1500 m), and ISO (1500–2000 m) patterns. Figure 7c displays the parallel (blue segment) and perpendicular (red segment) paths at the crossover. Each segment has 40 s of travel. The grid line data of these two segments were collected on the same day (12 August 2008) without any radar setting changes. The data were SAR (synthetic aperture radar) processed first, and then the multiple power–depth profiles within about 120 m distance for each segment were averaged to obtain the stacked power–depth profiles that were normalized to the maximum power of the perpendicular segment. The cross mark in Fig. 2a shows the location of the crossover related to the four circles. Figure 7d and e present the details of the power–depth profiles of the ANISO-1 and ANISO-2 patterns, respectively. For the power–depth profiles of ANISO-1 pattern, the power received along the segment perpendicular to the ice divide is constantly higher than that along the segment parallel to the ice divide, and the power difference is 1.5 dB on average, with a maximum difference of 6.3 dB at the depth of $\sim \mathrm{462}$ m. For the power–depth profiles of ANISO-2 pattern, the power received along the segment perpendicular to the ice divide is constantly lower than that along the segment parallel to the ice divide, and the power difference is 1.9 dB on average, with a maximum difference of 4.4 dB at the depth of $\sim \mathrm{1459}$ m.

Comparison with ice core measurements

The COF along the full NEEM ice core was measured with thin vertical sections of the core in field by an automatic ice texture analyzer every 10 m from a depth of 33 m down to 2461 m (Weikusat and Kipfstuhl, 2010), presented using the second-order orientation tensor ${\mathit{a}}^{\left(\mathrm{2}\right)}$ (Eichler et al., 2013; Montagnat et al., 2014). The eigenvalues of ${\mathit{a}}^{\left(\mathrm{2}\right)}$ represent the lengths of the ellipsoid axes best fitting the density distribution of the ice crystal $c$ axis orientations; its eigenvectors give the directions of the axes of the ellipsoid. Figure 8 displays the eigenvalues of ${\mathit{a}}^{\left(\mathrm{2}\right)}$ side by side with the radar echogram of HV polarization measurements at Circle 3 adopted from Fig. 3b for a visual comparison (the radar echogram is detrended to reduce the dynamic range for clearly displaying the weak layers close to the ice bed). According to Fig. 8, the COF around NEEM can be characterized with five sections (A, B, C, D, and E) along the whole ice column from surface to the bed. For Section A (From $\sim$34 m down to $\sim \mathrm{400}$ m), the COF is the isotropic single maximum ($a_{\mathrm{3}}^{\left(\mathrm{2}\right)}\approx a_{\mathrm{2}}^{\left(\mathrm{2}\right)}<\mathrm{1}/\mathrm{3}$) and we cannot observe the HV power distinction pattern. The COF evolves towards an elongated single maximum ($a_{\mathrm{3}}^{\left(\mathrm{2}\right)}<a_{\mathrm{2}}^{\left(\mathrm{2}\right)}<\mathrm{1}/\mathrm{3}$) from shallow to deep depths in Section B with an anisotropic fabric for depths down to $\sim \mathrm{1400}$ m, and the HV power distinction pattern becomes more obvious as the depth increases. The sudden drop both in COF eigenvalues and radar signal intensity in Section C between depths of $\sim \mathrm{1400}$ to $\sim \mathrm{1500}$ m corresponds to the Wisconsin–Holocene climate transition and results in a strengthened and sharpened single maximum COF ($a_{\mathrm{3}}^{\left(\mathrm{2}\right)}\approx a_{\mathrm{2}}^{\left(\mathrm{2}\right)}\ll \mathrm{1}/\mathrm{3}$) in Section D between $\sim \mathrm{1500}$ and $\sim \mathrm{2200}$ m. The HV power distinction pattern of birefringence continues from Section B through C and D to the three distinct deep ice layers at 1753, 1832, and 1894 m. In Section E from $\sim \mathrm{2200}$ m to the bottom, the ice is stratigraphically folded. The COF is multi-clustered in certain layers. Beyond 1894 m, we observe weak reflectors close to the ice bed from radar measurements.

(a) Fabric profile along the NEEM ice core represented by the eigenvalues of the second-order orientation tensor (yellow: $a_{\mathrm{1}}^{\left(\mathrm{2}\right)}$, red: $a_{\mathrm{2}}^{\left(\mathrm{2}\right)}$, and blue: $a_{\mathrm{3}}^{\left(\mathrm{2}\right)})$; Eichler et al., 2013; Montagnat et al., 2014); (b) Radar echogram of HV measurements at Circle 3, adopted from Fig. 3b for visual comparison with the ice core COF. A: single max; B: elongated single max; C: transition from B to D; D: strengthened single max; E: folded ice.

[Figure omitted. See PDF]

The HH power variation patterns revealed by Fig. 7 match with the second-order orientation tensor ${\mathit{a}}^{\left(\mathrm{2}\right)}$ shown in Fig. 8a. Roughly the pattern is ISO in Section A, ANISO-1 in Section B, ANISO-2 in Section C and ISO in Section D. The orientation of the crystals in the core is only known relative to the vertical axis. The orientation of the core in the $x$–$y$ plane after it was extracted is not known. Therefore, direct checks of ANISO-1 and ANISO-2 are not possible with the core data. In the following paragraphs, we compare the vertical orientation information revealed by radar polarization measurements and ice core data.

Figure 9 includes the pole figures (stereographic projections) projected into the sample plane (vertical) for depths at 330, 561, 849.8, and 1025.8 m, in which the black squares represent the eigenvectors of ${\mathit{a}}^{\left(\mathrm{2}\right)}$ and the figure the associated eigenvalues. In the pole figures, the center of the circle corresponds to a colatitude of zero. The colatitude is $\mathrm{90}{}^{\circ }$ and the azimuth is between $\left[\mathrm{0}\mathrm{2}\mathit{\pi }\right]$ along the circle. In Fig. 9, we can visually observe that the majority of ice crystals cluster around the vertical axis with a small offset tilt angle.

From Eq. (4), the wave number difference between the ordinary and extraordinary waves in birefringent ice is $$k_{\mathrm{2}}-k_{\mathrm{1}}=\mathrm{\varnothing }/\left(\mathrm{2}z\right).$$ From Eq. (7) in Hargreaves (1977), the anisotropic dielectric permittivity ${\mathit{ϵ}}_{z_{\mathrm{p}}}-{\mathit{ϵ}}_{x_{\mathrm{p}}}$ is related to the wave number difference by $$\begin{array}{cc}{\displaystyle }& {\displaystyle }({\mathit{ϵ}}_{z_{\mathrm{p}}}-{\mathit{ϵ}}_{x_{\mathrm{p}}}){\sin }^{\mathrm{2}}\mathit{\beta }={\displaystyle \frac{c^{\mathrm{2}}}{{\mathit{\omega }}^{\mathrm{2}}}}\left(k_{\mathrm{2}}^{\mathrm{2}}-k_{\mathrm{1}}^{\mathrm{2}}\right)\\ & {\displaystyle }& {\displaystyle }\approx (k_{\mathrm{2}}-k_{\mathrm{1}})(\mathrm{2}n/k_{\mathrm{0}}),\end{array}$$ where $\mathit{\beta }$ is the angle between the principal axis $z_{\mathrm{p}}$ and the $z$ axis, $n=\mathrm{1.77}$ is the ice refractive index, and $k_{\mathrm{0}}=\mathit{\pi }$ m${}^{-\mathrm{1}}$ is the wave number in air at 150 MHz. By substituting Eq. (9) into Eq. (10) we get $$\mathit{\beta }={\sin }^{-\mathrm{1}}\left[\sqrt{{\displaystyle \frac{n\mathrm{\varnothing }}{z{k}_{\mathrm{0}}\left({\mathit{ϵ}}_{z_{\mathrm{p}}}-{\mathit{ϵ}}_{x_{\mathrm{p}}}\right)}}}\right].$$

Fabric represented in stereographic projection/pole figure, projected into the vertical planes of the thin sections at depths of 330, 561, 849.8, and 1025.8 m. The OY direction in the pole figure is approximately the vertical ice core axis. The center of the diagram corresponds to the approximate horizontal. Only 10 000 data points randomly selected are plotted for each thin section. The color bars correspond to the density of plotted data points.

[Figure omitted. See PDF]

For small tilt angle $\mathit{\beta }$, we can obtain an approximate biaxial model described by Eq. (4) in Fujita et al. (2006) by rotating the uniaxial tensor about the principal axis $x_{\mathrm{p}}$ by $\mathit{\beta }$.

Using the inferred birefringence phase shift from the radar multi-polarization measurements shown in Fig. 6a, for ${\mathit{ϵ}}_{z_{\mathrm{p}}}-{\mathit{ϵ}}_{x_{\mathrm{p}}}=\mathrm{0.034}$, the average angle of the optic axis with respect to the vertical axis at different depths is estimated according to Eq. (11) and the result is presented in Fig. 10 at the locations of Circle 1 to Circle 4. The mean value of all the points in Fig. 10 is $\mathrm{11.6}{}^{\circ }$ with a standard deviation of $\mathrm{0.8}{}^{\circ }$, and the minimum and maximum are about $\mathrm{10}$ and $\mathrm{14}{}^{\circ }$, respectively.

By applying the concept of the effective medium (Maurel et al., 2016), we calculate the effective colatitude at 24 depths between 330 and 1025.8 m using the azimuth and colatitude measurements of the $c$ axis of ice crystals at these depths. The orientation of the $c$ axis for a single crystal is represented by $$c=(c_{\mathrm{1}},c_{\mathrm{2}},c_{\mathrm{3}})=(\sin \mathit{\theta }\cos \mathit{\phi },\sin \mathit{\theta }\sin \mathit{\phi },\cos \mathit{\theta }),$$ where $\mathit{\phi }\in \left[\mathrm{0}\mathrm{2}\mathit{\pi }\right]$ and $\mathit{\theta }\in \left[\mathrm{0}\frac{\mathit{\pi }}{\mathrm{2}}\right]$ are the azimuth and colatitude, respectively. The third component of the effective $c$ axis can be calculated by $$c_{\mathrm{3}}^{\mathrm{eff}}=\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{\mathrm{0}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}{p}_{\mathrm{n}}p(\mathit{\phi },\mathit{\theta })\cos \mathit{\theta }\sin \mathit{\theta }\mathrm{d}\mathit{\theta }\mathrm{d}\mathit{\phi },$$ where $p_{\mathrm{n}}$ is a normalization factor and $p(\mathit{\phi },\mathit{\theta })$ is the $c$ axis probability distribution function to meet the following normalization condition: $$\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{\mathrm{0}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}{p}_{\mathrm{n}}p(\mathit{\phi },\mathit{\theta })\sin \mathit{\theta }\mathrm{d}\mathit{\theta }\mathrm{d}\mathit{\phi }=\mathrm{1}.$$ Combining Eqs. (13) and (14) we get $$c_{\mathrm{3}}^{\mathrm{eff}}={\displaystyle \frac{{\int }_{\mathrm{0}}^{\mathrm{2}\mathit{\pi }}{\int }_{\mathrm{0}}^{\mathit{\pi }/\mathrm{2}}p(\mathit{\phi },\mathit{\theta })\cos \mathit{\theta }\sin \mathit{\theta }\mathrm{d}\mathit{\theta }\mathrm{d}\mathit{\phi }}{{\int }_{\mathrm{0}}^{\mathrm{2}\mathit{\pi }}{\int }_{\mathrm{0}}^{\mathit{\pi }/\mathrm{2}}p(\mathit{\phi },\mathit{\theta })\sin \mathit{\theta }\mathrm{d}\mathit{\theta }\mathrm{d}\mathit{\phi }}}.$$ We also get $${\mathit{\theta }}^{\mathrm{eff}}={\cos }^{-\mathrm{1}}\left(c_{\mathrm{3}}^{\mathrm{eff}}\right).$$

Comparison between the average angles of optic axis estimated from multi-polarization measurements and the effective $c$ axis angles estimated from ice core data at different depths.

[Figure omitted. See PDF]

Figure 11 displays an example of the measured $c$ axis orientation probability distribution as a function of $\mathit{\phi }$ and $\mathit{\theta }$ at a depth of 330 m. The calculated effective colatitudes are listed in Table 2, derived by numerically integrating Eq. (15) with these measurements at the 24 depths (during the numerical integration, a step of $\mathrm{1}{}^{\circ }$ was used for both d$\mathit{\phi }$ and d$\mathit{\theta }$). The average of these effective colatitudes is $\mathrm{80.4}{}^{\circ }$ with a standard deviation of $\mathrm{1.2}{}^{\circ }$. This corresponds to an effective $c$ axis tilt angle of $\mathrm{9.6}{}^{\circ }$ from the vertical axis, in very good agreement with the averaged angle $\mathit{\beta }=\mathrm{11.6}{}^{\circ }$ between the optic axis and the vertical axis derived from the radar multi-polarization measurements. The effective colatitudes at the 24 depths are converted to the tilt angle from the vertical and overlapped in Fig. 10 (magenta circles) for a straightforward visual comparison with the previously estimated ice optic axis tilt angle $\mathit{\beta }$. The minimum and maximum effective $c$ axis tilt angles from the vertical axis are $\sim \mathrm{5.5}$ and $\sim \mathrm{11.7}{}^{\circ }$ at 561 and 849.8 m, respectively.

The effective colatitude may not be physically and completely equivalent to the angle between the ice optic axis and the vertical axis, but the agreement shown between them may provide some insights on how ice COF will affect the radio wave propagations in the medium. The orientation of the thin sections relative to the main direction of the ice core was controlled with an accuracy of about $\mathrm{1}{}^{\circ }$; the orientation of the ice core by itself (relative to the vertical) is never well known and could vary by as much as $\mathrm{3}{}^{\circ }$. These experimental errors and the limited number of sampled crystals may affect the agreement between the averaged $\mathit{\beta }$ and ${\mathit{\theta }}^{\mathrm{eff}}$.

Measured $c$ axis orientation distributions as a function of $\mathit{\phi }$ and $\mathit{\theta }$ at depth 330 m. Any measurements with a quality factor of less than 70 were disregarded and not included in the distribution plots (see Peternell et al., 2010, for details about the quality factor). In total, there are qualified 2 032 533 crystal $c$ axis orientation measurements at this depth.

[Figure omitted. See PDF]

Effective colatitudes ${\mathit{\theta }}^{\mathrm{eff}}$ (${}^{\circ }$) calculated from COF measurements at different depths

Depth (m)	330	352	385	396	418	440	462	506	528	561	583	607.3
${\mathit{\theta }}^{\mathrm{eff}}$ (${}^{\circ }$)	80.4	81.7	79.9	80.4	79.9	81.1	79.8	80.2	80.2	84.5	79.7	80.9
Depth (m)	652.4	674.4	717.8	739.8	755.8	783.8	822.3	849.8	888.8	921.8	998.2	1025.8
${\mathit{\theta }}^{\mathrm{eff}}$ (${}^{\circ }$)	80.8	81.3	80.2	79.1	80.9	80.5	79.4	78.3	80.8	79.8	83.0	80.3

Summary, conclusions, and discussion

In our multi-channel and polarimetric radar measurements performed around the NEEM ice core site in 2008, we identified a power reduction of about 12.7 dB, transmitted towards the nadir in V-polarization because of the effects of the structure truss, mutual coupling on the antenna radiation pattern, and other unknown reasons. Despite this issue, redundant HH and HV polarization measurements were successfully completed, providing a valuable data set that reveals ice birefringence and COF characteristics. We presented a systematic method to process the data, effective for reducing speckle artifacts, enhancing the signal-to-noise ratio of ice layers, reducing power fading, and studying the bulk ice sheet birefringence and COF as a function of depth. The major findings from our analysis are as follows.

• The comparison of the power variation patterns for HH and HV measurements with simulated theory-based patterns suggests that the ice COF at NEEM site is an elongated single maximum for depths between 450 and 1500 m, which agrees with the COF measurements from the ice core.

• The birefringence phase shifts between the ordinary and extraordinary waves increase nonlinearly as a function of depth.

• The ice optic axis lies in vertical planes that are approximately perpendicular to or align with the ice divide, depending on depth. The spatial variations of specific optic orientations are about $\pm \mathrm{10}{}^{\circ }$ across the 100 km${}^{\mathrm{2}}$ area around the NEEM site (see Fig. 12b). In general, the ice optic axis plane rotates slightly counterclockwise as depth increases.

• For the depth range between 324 and 1000 m, the ice optic axis has an average tilt angle of $\sim \mathrm{11.6}{}^{\circ }$ from the vertical axis, with a standard deviation of $\mathrm{0.8}{}^{\circ }$ (see Fig. 12c). The average effective $c$ axis tilt angle from the vertical axis for 24 depths between 330 and 1025.8 is $\sim \mathrm{9.6}{}^{\circ }$ with a standard deviation of $\mathrm{1.2}{}^{\circ }$. The agreement between the ice optic axis angle and effective $c$ axis tilt angle provides insight into how ice COF will affect the radio wave propagations in the medium.

• The anisotropic patterns and ice optic axis orientations revealed by the radar polarimetric measurements provide insight into the ice flow history of the NEEM site.

• The very good agreement between the radar observations, simulations, and ice core fabric data provides assurance of the effectiveness of polarimetric radio echo sounding techniques to infer profiles of ice fabric in locations where no ice core data are available.

Illustration of the major findings from the multi-polarization radar measurements. (a) Diverging ice flow scenario at ice divide; (b) sketch of girdle COF formed by diverging ice flow in which the point $z_{\mathrm{p}}$ is the projection of the ice optic axis, and rotation arrow indicates the spatial variations of the ice fabric orientation around the NEEM site; (c) the average tilt angle of ice optic axis $z_{\mathrm{p}}$ at the NEEM site.

[Figure omitted. See PDF]

The orientation of the ice COF is dictated by stresses in the ice column, and the preferred orientation of the girdle COF tends to be perpendicular to the direction of tension (van der Veen and Whillans, 1994). The results of this investigation show two different girdle orientations within the area of study. In the case of ANISO-1 (from $\sim \mathrm{450}$ to $\sim \mathrm{1065}$ m), the girdle COF is oriented perpendicular to the ice divide, indicating parallel ice flows along the ice divide (see Fig. 12a, the main ice flow velocity vector is parallel to the divide). Around the NEEM site, the main slope is in the divide direction, the slope is very low on both sides of the divide, and thus the flow scenario is divergent from the main flow along the divide (Fabien Gillet-Chaulet, personal communication, 2018). In the case of ANISO-2 (from $\sim \mathrm{1375}$ to 1575 m), the girdle COF is oriented parallel to the ice divide, suggesting the ice once flowed away from the ice divide on both sides in a direction orthogonal to the ice divide. The presence of both COF orientations shows the complex flow history of the area of study and demonstrates the need for additional investigations in both present and previous ice flows of the region.

The multi-polarization radar measurements after pulse compression and multi-channel combining are available at https://data.cresis.ku.edu/data/rds/2008_Greenland_Ground/CSARP_qlook_polarimetric_chan_combined/ (last access: 16 April 2018).

The Supplement related to this article is available online at https://doi.org/10.5194/tc-12-2689-2018-supplement.

All authors contributed to this work. JL processed the radar data, performed the analysis of both radar and ice core records, and led the writing of the manuscript. JAVG performed the initial analysis of the radar data. CL contributed to the radar system design and collected the radar data. AH designed and developed the antenna array for the radar system. PG initiated the research project as the CReSIS principal investigator and contributed to the radar system design. MM and IW contributed to the ice core data collection and analysis. FRM contributed to the radar system design and implementation. JP contributed to the processing toolbox for the radar data. All co-authors provided input to the manuscript.

The authors declare that they have no conflict of interest.

Acknowledgements

This project was supported by NSF grant ANT-0424589 and NASA grant NNX16AH54G. The authors would like to thank all the students and staff at CReSIS who took part in the development of the radar system and the data collection in the field. They also acknowledge the support from the NEEM team during the data collection and particularly Sepp Kipfstuhl for conducting most of the COF measurements. Ilka Weikusat acknowledges funding by HGF VH-NG-802. Rachel McCarthy James is acknowledged for editing the manuscript. The authors would also like to thank the reviewers and editors for their valuable input to improve the paper.Edited by: Martin Schneebeli Reviewed by: Howard Conway and Zoe Courville