Summit of the East Antarctic Ice Sheet underlain by thick ice-crystal fabric layers linked to glacial–interglacial environmental change

Abstract Ice cores in Antarctica and Greenland reveal ice-crystal fabrics that can be softer under simple shear compared with isotropic ice. Owing to the sparseness of ice cores in regions away from the ice divide, we currently lack information about the spatial distribution of ice fabrics and its association with ice flow. Radio-wave reflections are influenced by ice-crystal alignments, allowing them to be tracked provided reflections are recorded simultaneously in orthogonal orientations (polarimetric measurements). Here, we image spatial variations in the thickness and extent of ice fabric across Dome A in East Antarctica, by interpreting polarimetric radar data. We identify four prominent fabric units, each several hundred metres thick, extending over hundreds of square kilometres. By tracing internal ice-sheet layering to the Vostok ice core, we are able to determine the approximate depth–age profile at Dome A. The fabric units correlate with glacial–interglacial cycles, most noticeably revealing crystal alignment contrasts between the Eemian and the glacial episodes before and after. The anisotropy within these fabric layers has a spatial pattern determined by ice flow over subglacial topography.

Ice sheets contain several types of single ice crystals, such as a perfect single pole, an elongated single pole and a vertical girdle crystal (Fig. 1). Ice crystals experience an anisotropic response to stress, whereby their principal axes, termed C-axes, rotate during deformation towards the direction of compression and away from the tensional direction (Fig. 2;Alley 1988;Azuma et al. 1999;Durand et al. 2009;Montagnat et al. 2014). Collectively, icecrystal alignment leads to macro-crystallographic orientations, or ice-sheet fabric, which is the product of both ice rheology and stress. As the strain rate in anisotropic ice under simple shear (Azuma 1994;Durand et al. 2007;Kennedy et al. 2013) can be up to an order of magnitude greater than in isotropic ice (Azuma 1994), a feedback is possible whereby ice fabrics lead to preferential deformation and, in turn, further anisotropy (Alley 1988;Azuma 1994). The spatial variation of ice rheology is fundamental to ice-sheet flow modelling, yet it is normally averaged through depth as an isotropic medium, as in Glen's flow 'law' (Glen 1954). Differences in the velocity of cold-based ice between observations (Thorsteinsson et al. 1999) and modelling (Gillet-Chaulet et al. 2006;Martin & Gudmundsson 2012), in which respective ice-surface elevations can vary by hundreds of metres, and ice velocities by tens of metres per year, may be attributable to the rheological impacts of ice fabric development. Crystal alignment in ice cores is analysed routinely, revealing marked and step-wise changes in the fabric of ice deposited across glacial-interglacial periods (Durand et al. 2007;Montagnat et al. 2014). Additionally, sonic logging (Gusmeroli et al. 2012) shows how variations in ice-sheet fabric is correlated highly with oxygen isotope ratios, implying a link with atmospheric parameters (Paterson 1991;Kennedy et al. 2013). The spatial distribution of ice-sheet fabric is largely unknown, however, as ice cores are heavily restricted in number and distribution across Antarctica and Greenland. As a consequence, the crystallographic structure of ice in large ice sheets, including its development through time and space, and its potential rheological consequences, remains uncertain.
Dome A is the highest ice divide on the East Antarctic Ice Sheet , with surface elevations reaching over 4100 m above sea-level. During the 2007-09 period of the International Polar Year the Dome A summit region was surveyed by an ice-penetrating polarimetric radar (Sun et al. 2009;Fig. 3). Initial inspection of the non-polarized radar data revealed the first detailed images of the Gamburtsev Subglacial Mountains, which were found to contain wellpreserved upland glacial geomorphological features and a clearly defined water drainage pattern originating from fluvial activity in preglacial times (Sun et al. 2009). Numerical ice-flow modelling suggests these mountains were the nucleation site of the East Antarctic ice sheet, and are likely therefore to have been continuously glaciated for at least 14 Ma (DeConto & Pollard 2003).
Because radio-waves reflect off ice fabrics as a function of orientation, polarimetric radars, which measure radio-wave reflections simultaneously in orthogonal axes, are able to detect gross levels of ice-crystal orientations. Here, we analyse polarimetric data from the Dome A radar survey, revealing insights into the deep-ice crystal structure at and around the ice divide, and its relation with subglacial topography. By linking the data stratigraphically to the Vostok ice core, using airborne radar data, we also provide an assessment of the temporal evolution of ice fabric.

Typical fabrics in ice sheets
Deep ice cores, and some laboratory measurements confirm a direct association between ice-fabric changes and their permittivity (Matsuoka et al. 1997;Fujita et al. 2000Fujita et al. , 2006Drews et al. 2012). Figure 1 and Table 1 show the usual fabric types and their permittivity tensors. Electric permittivity is a measure of how an electric field affects, and is affected by, a dielectric medium. In tri-dimensional anisotropic media, the permittivity is a tensor relating the three components of the electric displacement with those of the electric field. The permittivity tensor can be expressed as a diagonal tensor, as follows: (1) We assume here that ε 1 ≤ ε 2 ≤ ε 3 . In the case of a random ice-crystal fabric, the response of the electric field within the medium is independent of its orientation and ε 1 = ε 2 = ε 3 . If ε 1 = ε 2 = ε ⊥ and ε 3 = ε // , an idealized perfect singlepole behaviour occurs (Fig. 1b). For a single-pole ice-fabric distribution, as the ice crystals are aligned in planes, distinct electrical properties will exist along and across those planes. If we denote ε ⊥ and   ε // the permittivity perpendicular and parallel to the normal vector of the crystal planes, respectively, the permittivity tensor can be expressed as ε 1 = ε 2 = ε ⊥ and ε 3 = ε // (Fig. 1b).
If the ice crystals have elongated single poles then ε 1 < ε 2 < ε 3 (Fig. 1c). The ellipse ratio in this fabric corresponds to the anisotropy of permittivity in the horizontal direction. The larger this ratio is, the higher the level of anisotropic permittivity.

Polarimetric radio waves and Stokes parameters
Polarimetric radars have orthogonal transmitting and receiving antennas, receiving four series of polarized echo (noted as HH, HV, VH and VV). Relative echo powers from these configurations vary through anisotropic media. The classic Stokes vector parameters (referred to below as S 0 , S 1 , S 2 and S 3 ) can be used to evaluate and quantify the polarization of such radio-waves as follows.
The Stokes vector parameters, defined by George Gabriel Stokes in 1852, are widely used in electromagnetics and optics for describing the polarization state of electromagnetic radiation, such as its total intensity and degree of polarization and the shape parameters of the polarization ellipse (Stokes 1852;Chandrasekhar 1960;Doake 1981;Mott 2006). As H and V are orthogonal, we can regard them as axes to build the corresponding Stokes parameters, which can be represented as a vector as follows: Here where E H and E V are the components of electrical field from the orthogonal receiver antennae; φ is the phase difference between E H and E V ; and Re and Im are the real and imaginary parts of an electrical field E. The normalized Stokes vector parameters can effectively eliminate the influence of transmitting power, antenna interval and so on, and is represented as follows: While an electromagnetic wave propagates in an anisotropic media, it can induce both ordinary waves and extraordinary waves in orthogonal directions, each with different propagating velocities. Therefore the incident wave in the H or V direction can induce reflections both in the same direction and orthogonally. We can describe the relationship of reflection and incident waves as: Here, E represent the electrical field vector. Its superscript 'i' means that it is the incident wave and 'r' indicates the reflection wave. R is a 2 × 2 reflection coefficient matrix describing the relationship between the incident and reflecting waves.  Figure 1 Type of fabric The main components of the permittivity tensor are shown in both the primary axis and in a Cartesian system of coordinates (x, y, z), where z represents the vertical.
For induced reflections recorded in the H antenna, we can obtain: Substituting equation (9) into equation (2), we get the following four relations: For inducing reflections measured in the V antenna we obtain similar expressions.
S 0 refers to the total echo-power from the orthogonal antennae, and resembles the reflected power using in the conventional sole-channel radar system. S 1 , S 2 and S 3 provide the difference in reflected power between the orthogonal electric fields in the standard Cartesian arrangement, a Cartesian system rotated by 45°and in a circular coordinate arrangement, respectively, thus providing measures of the anisotropy in ice permittivity.
Several authors have confirmed the relationship between the echo difference in orthogonal directions and the permittivity tensor of ice (Hargreaves 1977;Fujita et al. 2006;Matsuoka et al. 2012), which is linked to the second-order tensor of ice fabric (Matsuoka et al. 2009). Hence, in the absence of plausible alternative explanations, it can be assumed that the Stokes vector parameter description of the polarimetric radar reflections S 1 , S 2 and S 3 expresses ice-sheet anisotropy caused by ice crystal fabric.

Radar survey
The polarimetric radar used a frequency of 179 MHz and recorded radio waves in orthogonal orientations over a square grid of fourteen 30 km-long transects with along-track sample acquisition every c. 17 m (Fig. 3). The radar employed dual receivers in orthogonal orientations (H and V ) with a timesharing sole transmitter (H or V ). It comprised four assembly antenna pairs: HH, HV, VH and VV. An in-phase/quadrature (I/Q) demodulator can measure the amplitude, phase and instantaneous frequency of radar signals, so it is widely used in polarimetric radar systems. Normally each antenna pair records both I and Q components for detecting phase. Therefore there are eight data records in the polarimetric radar files necessary to measure information of both the amplitude and phase of reflected waves. Further details of the radar are provided in Sun et al. (2009).
The I and Q components of the four antenna pairs were imported from raw data files, from which the Stokes vectors of both H and V transmitters were calculated. The four parameters in the Stokes vector each have an explicit physical meaning. S 0 represents the average power of the EM waves. S 1 provides the echo difference between H and V components. Its value can be positive, zero or negative depending on whether linear polarization of H or V is dominating. S 2 provides the echo difference of orthogonal components between 45 and −45°. S 3 represents the echo difference of orthogonal components within a circular coordinate system. After normalization, S 1 , S 2 and S 3 all range between −1 and +1.
S 0 has similar results to traditional (nonpolarimetric) radar reflections. It can display clear isochrones and subglacial topography. In contrast S 1 , S 2 and S 3 display obvious bands within the normal isochrone arrangements (e.g. red (positive) and blue (negative) layers in Fig. 4). The units in S 1 and S 2 are more distinct than in S 3 and show little difference between the H and V transmitters. By adding S 0 , S 1 and the average power from H and V together, the positive layers link to reflections recorded in the V component, and the negative layers link to the H component. If the coordinates are fixed (i.e. as in a straight survey line), variation of the Stokes vectors in the H and V components within different layers of ice must be attributed to differences in the ice crystal fabric (i.e. different layers have different fabric and orientation tensors).
The Stokes vectors were calculated in each of the transects (Fig. 4), revealing four obvious englacial layered units in which discrete variations in S 1 were noticed, named T 1 , T 2 , T 3 and T 4 , respectively (T 0 designates the bed). Beneath T 4 , further units of alternating S 1 occur, but they are less obvious and not as easy to trace owing to the lack of returned power from deeper englacial reflectors. Broadly, positive S 1 units (red-yellow zones in Fig. 4) correspond to high levels of anisotropy in the horizontal plane (as the radar is insensitive to the anisotropy in the direction of transmission), and negative S 1 units (blue zones) are derived from either isotropic ice or low levels of anisotropy in the horizontal plane. As expected, the deeper S 1 ice fabric units thin gradually with depth owing to vertical compression (Fig. 4). Unit T 1 is the thickest of those measured, from the surface to ice depths of c. 800-1000 m. Based on crystallographic information from deep ice cores, this upper unit will have a nearly random (isotropic) fabric (Azuma 1994;Azuma et al. 1999;Fujita et al. 2006;Durand et al. 2009;Montagnat et al. 2014). Unit T 2 (positive S 1 ) is up to c. 400 m thick and contains an abundance of highly anisotropic ice related to a vertical girdle fabric distribution. T 3 (negative S 1 ) is a product of low levels of horizontal anisotropy, as the fabric clusters towards a broad single-pole distribution (Hargreaves 1977;Eisen et al. 2007;Matsuoka et al. 2009). T 4 (highly positive S 1 ) is similar in crystal arrangement to T 2 . While both the thickness and depth of the units fluctuate horizontally substantially, we estimate that c. 20% of the ice column beneath Dome A, and up to c. 50% of the ice column at depths >1 km around this ice dome, are characterized by high levels of anisotropy (primarily from T 2 , T 4 and T 6 ; Fig. 4).

Link to Vostok ice core
We are able to link radar isochrones from Kunlun Station to the Vostok ice core, allowing us to estimate a depth-age relationship at Dome A (Figs 5 & 6). The depth of the ice fabric units at Kunlun is then used to assign estimated ages to the fabric units (Table 2). We also used the stratigraphic connection to Vostok to translate the ice core dust and ΔT records to Dome A (Fig. 7).
Isochrone layers were picked along two airborne radar transects, collected by the Alfred Wegener Institute. The radar used operates at 150 MHz, alternating between 600 and 60 ns bursts, with internal stacking of shots prior to storage and a pulse repetition frequency of 20 Hz. The effective distance between 60 ns shots is c. 8 m, and for interpretation the shots are stacked two-fold (so that c. 8 m becomes c. 16 m). The two-way travel times of the layers were converted into ice depth using the Kovacs formula (Kovacs et al. 1995) and ice density information (Barkov 1970;Ekaykin et al. 2010), with an accuracy of c. 5%. The radar data were differentiated and amplitude balanced to aid the extraction of the layers. The age-depth relationship at Vostok Station was taken from Parrenin et al. (2004). The layers were traced using a commercial software package from Landmark, which displays picked horizons at cross-over points and thus allows an evaluation the quality of the data to ensure the correct isochrone is identified along the transect. Importantly, layers were traced along both transects from Dome A and compared at their intersection at Vostok Station. Agreement between the layers at Vostok Station provides confidence that the correct layers were followed to Dome A, and that the depth-age relationship there is reliable. The error for picking both the ice-surface reflection and an isochrone is about five sample points, which corresponds to 66.65 ns or c. 5.6 m in ice thickness.

Ages of fabric layers
From the stratigraphic correlation between Vostok and Dome A, we estimate that the base of T 2 is c. 105 kyr, which roughly corresponds with the end of the Eemian interglacial period. T 3 is aged between c. 105 and c. 136 kyr, and there is a good correspondence between the peak of S 1 within T 3 and the peak of the Eemian. T 4 is between c. 136 and c. 190 kyr, similar to MIS 6. Two faint units below T 4 , named T 5 and T 6 , are possibly consistent with previous glacial/interglacial periods.
Broadly, positive S 1 units (T 2 , T 4 ) correspond well to glacial periods and negative S 1 units (T 3 , T 5 ) are associated with interglacials (Fig. 7). This suggests that gross crystal fabric across Dome A is influenced by past environmental change at the glacial-interglacial transition.

Crystal fabrics at Dome A: distribution and significance
In a vertical column of ice near the ice divide, the orientation of ice crystals will gradually tend to a predominant (vertical) orientation under compression. In doing so, the ice fabric develops from 'random' to 'single pole'. Horizontally, perfect single pole and 'vertical girdle' crystal alignments represent the two end-member extreme fabrics (Figs 1 & 2a) that could be expected in idealized ice-divide and ice-sheet margin regions, respectively. The resulting development of ice fabric in layers T 1 -T 5 is illustrated in Figure 8. An alternative, theoretical, explanation is that under shear the crystal axis experiences rotation in line with the stress field (Fig. 2b). Both Fig. 6. Depth-age curve at the Dome A summit. The blue points denote the age-depth data from dated radar internal horizons tracked from the Vostok ice core site (Fig. 5). The line is the 'best-fit' (polynomial) curve to the data points, and serves as an approximate depth-age relationship at Dome A.  Table 3, would yield similar polarimetric radar measurements. It is known that soft layers (glacial ice) are more easily deformed under stress than hard layers (interglacial ice) and isotropic ice (e.g. Paterson 1991; Thorsteinsson et al. 1999;Durand et al. 2007;Matsuoka et al. 2012;Kennedy et al. 2013). The deformation leads to the glacial ice layers more easily experiencing crystal realignment (Fig. 2a), which is what is recorded by the Stokes vectors measurements (Fujita et al. 2006;Wang et al. 2008;Matsuoka et al. 2012). Hence, the polarimetric radar Fig. 7. Comparison between ice fabric units and the Vostok dust and ΔT records. We construct a new age-depth relation at Dome A by dating radar internal horizons by tracking them from Kunlun Station to the Vostok ice core site (Fig. 6). T 2 -T 4 refer to the same fabric units as presented in Figure 4d. T 5 and T 6 are deeper units, and because both are faint we can only infer approximate ages. measures the thickness and distribution of layers with discrete ice fabrics.
The S 1 units align well among the individual transverses and cross profiles. We picked the top and bottom surfaces of all the S 1 units in each line (where S 1 became 0) to make three-dimensional wall and surface diagrams (Fig. 9), which display the extent and thickness of each fabric horizon revealing two obvious features. The first is significant S 1 variations in both ice depth and unit thickness. The second is changes to ice fabric caused by flow over macro-scale bed topography (Robin & Millar 1982). The fabric units are dominant across the left side of a major subglacial valley. On the right side of the subglacial valley, a cluster of relatively small peaks exists. Unit T 2 pinches out immediately above these peaks (Fig. 9). In this area, an arcuate belt with no S 1 signal covers the cluster of peaks. We infer that, above these smaller peaks, a pure compressive stress state (uniaxial or pure shear) dominates, causing the fabric to be characterized by a single-pole arrangement (Azuma 1994;Matsuoka et al. 2012). This induces isotropy in the horizontal direction, causing the S 1 to be near zero (c. 25 km mark in Fig. 4d), leading to local disappearance of the fabric units. Away from the peaks, simple shear will increase gradually, leading to the redevelopment of fabric units that can be seen on both sides (Fig. 9d).
Our polarimetric radar data reveal widespread ice-sheet anisotropy at Dome A. However, polarimetric radar only captures relative contrasts in electric permittivity between two orthogonal directions in the horizontal plane rather than the actual values. Interpretation of our radar data has been guided by typical ice-fabric distributions sampled in ice cores sites (e.g. Durand et al. 2007 at Dome C; Drews et al. 2012 at EPICA DML;and Montagnat et al. 2014 at NEEM). Ice-core data only show subtle changes in fabric during climatic transitions, however, in contrast with the marked transitions we image in radar data. This suggests that widespread but locally subtle transitions in ice fabric are likely to have a large effect on the permittivity of the ice sheet. Hence subtle transitions in ice fabric could potentially have a larger effect on ice flow dynamics than currently considered.
Our observations at Dome A are likely to hold for other ice divides, including central Greenland, as noted above. At Dome F in East Antarctica, using a dual frequency radar, simple-shear-dominated ice flow over large subglacial obstacles has been shown to lead to preferential recrystallization along bands of relatively soft ice (Fujita et al. 1999), meaning that the rheological consequences of initial anisotropy will evolve downstream. At Dome C, Siegert & Fujita (2001) analysed 60 MHz radar, and discovered that the form of the depth-to-power relationship of the returned radar signal could be used as a proxy for the along-track development of ice-crystal fabric. A similar investigation revealed how ice fabric may evolve between Ridge B and Lake Vostok (Siegert & Kwok 2000).
One situation where the downstream modification to ice crystal fabric and rheology may be highly relevant is in the preservation of major blocks of basal ice, interpreted to be caused by widespread and thick basal accretion (Bell et al. 2011). Stressinduced anisotropy-may lead to preferential englacial flow along bands, leading to smooth englacial flow over rough topography and, potentially, a partial flow 'decoupling' of the bulk of the ice sheet from basal ice units frozen to the bed. This new and alternative explanation for the so-called 'accretion ice' identified atop the Gamburtsev Subglacial Mountains (and other 'fuzzy' lower-englacial features seen in both Greenland and Antarctica) is consistent with the lack of significant surface uplift and velocity anomalies that a high rate of basal freezing (of the order of the surface rate of ice accumulation) should lead to, and is consistent with meteoric radar layering above the basal ice that drapes smoothly over it (see Wrona et al. 2017). Such flow may preserve the basal ice units as ancient ice blocks, which means that the rate of ice accumulation by freezing will be lower, and the age of ice considerably older (possibly by millions of years), than first thought.

Conclusions and further work
Polarimetric radar measurements, linked stratigraphically to the Vostok ice core, demonstrate that icesheet crystal alignment at the centre of large ice sheets is a function of both climate-induced impurities acting to 'soften' the ice during periods of glaciation and internal stresses within the ice sheet, acting to preferentially deform the soft ice, yielding ice-sheet anisotropy. These observations at and around Dome Table 3. Proposed dominant crystal fabrics at Dome A, following two possible crystal evolutions illustrated in Figure 2 Layer Option a in Figure 2 Option b in Figure 2 T   Table 3).
A are likely to hold for other ice divides, including central Greenland, where soft anisotropic ice has been measured in the lower portion of the NEEM ice core (NEEM Community Members 2013). The next steps to this area of research are twofold: first, to study ice crystal fabric developments using polarimetric radar around other ice divides, to understand the similarities to, and differences from, the situation at Dome A; and second, to investigate the degree to which the ice fabric uncovered at Dome A relates to modification in ice rheology. Icesheet modelling is capable of such analysis, through assimilation of three-dimensional internal layer structures into the ice flow parameter. This coupling between data and modelling is possible, but requires high-resolution radar data acquisition across study regions.