Polarimetric radar reveals the spatial distribution of ice fabric at domes in East Antarctica

domes in East Antarctica M. Reza Ershadi1, Reinhard Drews1, Carlos Martín2, Olaf Eisen3,5, Catherine Ritz4, Hugh Corr2, Julia Christmann3,6, Ole Zeising3,5, Angelika Humbert3,5, and Robert Mulvaney2 1Department of Geosciences, University of Tübingen, Tübingen, Germany 2British Antarctic Survey, Natural Environment Research Council, Cambridge, UK 3Alfred Wegener Institute Helmholtz-Centre for Polarand Marine Research, Bremerhaven, Germany 4University Grenoble Alpes, CNRS, IRD, IGE, Grenoble, France 5Department of Geosciences, University of Bremen, Bremen, Germany 6Institute of Applied Mechanics, University of Kaiserslautern, Germany Correspondence: Mohammadreza Ershadi (mohammadreza.ershadi@uni-tuebingen.de)

. The background color is the bed elevation (Morlighem et al., 2020). Yellow arrows are the magnitude and direction of the surface velocities for Dome C (Vittuari et al., 2004) and EDML (Wesche et al., 2007). The white strain ellipses mark the directions of the maximum and minimum strain rate. v1 and v2 are the ice-fabric's horizontal Eigenvectors, and they are based on the results in Sects. 4.1 and 4.2. Note that (a) and (b) have a different scale and orientation.

Quantitative metrics used to define the ice fabric
Ice crystallizes in the shape of hexagons and the direction normal to the basal plane is described with the c-axis (Hooke, 2005).
In a given strain regime, individual ice crystals orient themselves so that the bulk c-axis orientation of many crystals forms a 75 distinct pattern which we refer to as ice fabric. Elsewhere it is also described with Crystal Orientation Fabric (COF) or Lattice Preferred Orientation (LPO) (Weikusat et al., 2017). Development of ice fabric can lead to ice softening along with certain directions by up to a factor of 60 (Duval et al., 1983). The mechanical anisotropy is accompanied by a dielectric anisotropy to which the radio waves are sensitive described by the orientation tensor (Gödert, 2003;Gillet-Chaulet et al., 2006;Martín et al., 2009). The bulk ice-fabric pattern is often quantified using the Eigenvectors (v1,v2,v3) and Eigenvalues (λ1, λ2, λ3) of an 80 ellipsoid that best represents the c-axis orientation of all ice crystals in the sample. The Eigenvalues are normalized 6 https://doi.org/10.5194/tc-2020-370 Preprint. Discussion started: 20 January 2021 c Author(s) 2021. CC BY 4.0 License. and to be consistent with the past polarimetric radar studies, we assume λ1 < λ2 < λ3.

Data collection 90
The radar data in this study were collected using a phase-sensitive frequency-modulated continuous-wave radar system (Brennan et al., 2014;Nicholls et al., 2015) with a 200 MHz bandwidth and f c = 300 MHz center frequency. This radar emits linearly polarized electromagnetic waves using a slot antenna where the direction of the polarization plane is aligned with the direction of the electric field vector (e) in the antenna as shown in Fig. 2a.
We use terminology from satellite radar polarimetry to distinguish the directions of the polarization with H and V, although, 95 in a nadir-looking geometry, these are arbitrarily determined because H and V both have horizontal polarization plane at depth.
Here, we name the horizontal (H) and vertical (V) polarization plane consistent with Jordan et al. (2019). However, we want to point out that this definition is different to the one applicable to seismic shear waves, where vertically receiver (thus having a vertical component upon reflection at depth for non-vertical angles of incidence), and vice-versa. polarization planes are defined so that H is parallel to TR. v1 and v2 are the directions of the ice-fabric horizontal principal axes. θ is the angle between H and v1, and α is used for georeferencing. The model coordinate system is shown in Fig. 2c. The aerial line (TR) connects transmitter (Tx) and receiver (Rx), and by 100 convention, we assume that H is parallel to TR. v1, and v2, are the horizontal Eigenvectors which align with the direction of the smallest ( x ) and largest ( y ) horizontal principal permittivity, respectively (Fujita et al., 2006;Jordan et al., 2019). Hence, θ = 0°if H is aligned with v1. The angle α is measured by compass with ±15°uncertainty for georeferencing the data. Here, we use polar stereographic coordinates where anticlockwise rotation is positive.
Radar data at all the sites were collected at a fixed α, obtained from different antenna orientation in co-polarization (HH, VV) 105 and cross-polarization (HV, VH) configurations (Hargreaves, 1977;Fujita et al., 2006) as shown in Fig. 2b. We refer to these measurements as quad-polarimetric measurement. Radar data at Dome C were collected at 20 sites in January 2014. One of the sites is located within walking distance of the ice-core site EDC. The remaining 19 sites (termed E(ast)0-E18, and W(est)0.5-W18, with the numbers relating to the distance in km away from the dome) are aligned in a profile which is approximately perpendicular to the long axis of the dome and parallel to the flowline (Fig. 1a). At EDML, data were collected in January

Background of radar polarimetry
Radio signal propagation through ice sheets is polarization-dependent because of the dielectric anisotropy of the ice fabric.
If the direction of v3 is vertical and the remaining two Eigenvectors (v1, v2) are in the horizontal plane, then the relation 115 between the depth profile of the dielectric permittivity tensor and the orientation tensor is given by Fujita et al. (2006): For the dielectric permittivity at radio frequencies perpendicular to c-axes, we use ⊥ = 3.15 (Fujita et al., 2000), which is slightly lower than the value found by Bohleber et al. (2012). The value of a dielectric anisotropy for a single crystal is set to ∆ ≈ 0.034 (Matsuoka et al., 1997). The vertical v3 assumption in this study is justified through measurements at the EDC 120 ice core where the direction of v3 varies only by about 5°around the vertical (Durand et al., 2009). Elsewhere in ice sheets, this may not be the case, which will cause an additional source of horizontal birefringence (Matsuoka et al., 2009;Jordan et al., 2019).
We will model radio-wave propagation through birefringence ice using the method developed by Fujita et al. (2006). It includes transmission and reflection of initially linearly polarized waves emitted with two polarization modes (H and V, with 125 direction defined in the previous section). If z is the depth from the surface (positive downward), it assumes stratified ice with i = 1, ...N layers predicting the radar response as a function of the emitted polarization plane and ice-fabric parameters. Radar transmission (T) and reflection (Γ) are represented by 2 × 2 matrices only because radar signal propagation is insensitive to the vertically directed v3. The transmitted radar wave (e T ) and the corresponding radar reflection (e R ) are 2 × 1 vectors, with each component containing the electric field information of the H and V polarization components, respectively. Because only 130 relative phase and amplitude variations are considered, all information about the radio wave transmission and reflection can be inferred from the scattering matrix (S) at layer N : containing the complex scattering unit: where s HH and s V V are the co-polarized scattering signals, s HV and s V H the cross-polarized scattering signals, j is the imaginary unit, and k 0 = 2πf c c −1 0 is the wavenumber in vacuum with c 0 the speed of light in vacuum. A rotation between TR and v1 of the ice fabric in layer i, (θ i ), is accounted for by the rotation matrix R and its transpose (R ) Transmission of the signal is described by the transmission matrix T along the ice-fabric horizontal principal axes. T is a function of wavenumbers (k x , k y ), whereas the wavenumbers can be expressed as a function of dielectric permittivities ( x , y ) and electrical conductivities (σ x , σ y ) (Fujita et al., 2006).
145 where 0 and µ 0 are the dielectric permittivity in vacuum and the magnetic permeability in vacuum, respectively, and ω is the angular frequency. In this study we follow Fujita et al. (2006) and assume isotropic electrical conductivity (σ x = σ y ). Using Eq. (3), T can be written as a function of Eigenvalues: where it tracks the relative phase shifts induced by the dielectric anisotropy along the ice-fabric principal axes. The reflection 150 matrix Γ (Ulaby and Elachi, 1990) describes the reflection of the radio waves at an interface with changing dielectric properties containing the Fresnel reflection coefficients Γ x and Γ y calculated from depth-variable changes in permittivity (Paren, 1981;Ulaby and Elachi, 1990). In our case, these are caused by abrupt changes in ice-fabric orientation and/or magnitude. The 155 reflection ratio r = Γy Γx is a measure for the polarization dependence of the reflection boundary. If anisotropic scattering is caused by depth variable ice fabric only, then the reflection ratio at the interfaces i and i + 1 can be approximated (Paren, 1981;Drews et al., 2012) by: This offers at least in theory the option to fully reconstruct λ1, λ2, and λ3 in a nadir geometry, which will resolve the ice-fabric 160 types ambiguity as explained in Appendix B. Further details about the radar forward model implementation and definition of all the parameters in Eq. (6) are described in Fujita et al. (2006).
The parameters of interest that we aim to infer from the radar observations for each layer are the horizontal anisotropy ∆λ = λ2 − λ1, the ice-fabric orientation angle θ, and the reflection ratio r. All of these quantities may vary with depth. Much information is gained by interpreting the phase difference between the ordinary and extraordinary radio waves along v1 and 10 https://doi.org/10.5194/tc-2020-370 Preprint. Discussion started: 20 January 2021 c Author(s) 2021. CC BY 4.0 License.
where ∆φ x and ∆φ y are the phase shift at the reflection boundary, which we set both to zero (Fujita et al., 2006;Jordan et al., 2019Jordan et al., , 2020. The phase difference in Eq. (13) can vary rapidly with depth. Therefore, it is replaced by the coherence phase difference between s HH and s V V , which is a crucial development in the works from Dall (2010). The coherence phase 170 difference φ HHV V is the argument of the complex polarimetric coherence C HHV V , estimated via a discrete approximation, where M is the number of range bins used for vertical averaging, and b is the summation index. Using the depth gradient of φ HHV V in a Taylor expansion of Eq. (13) termed scaled phase derivative, which provides a way to relate the local phase gradient to ∆λ at the direction of the horizontal principal axes (Jordan et al., 2019(Jordan et al., , 2020 ∆λ(z) = Ψ(z, θ = 0°, 90°).
The ApRES stores the de-ramped signal (Brennan et al., 2014;Jordan et al., 2020), which is not represented in Eqs. (14) and

180
(15). The deramping corresponds to a complex conjugation of C HHV V (Jordan et al., 2020). Therefore, we use Eq. (14) for the models and the conjugate of Eq. (14) for the radar data to calculate the coherence phase.
We simplified Eq. (6) to a single layer case (Appendix C) showing that the polarity of Ψ can differentiate the direction of v1 and v2 (Appendix D). If the coherence phase is based on the received signal, v2 is in the direction of Ψ > 0 (TR v2), and v1 is in the direction of Ψ < 0 (i.e., TR v1). When using observations, the depth gradient calculation of φ HHV V is inherently 185 difficult because any differencing scheme amplifies noise (Chartrand, 2011). We follow Jordan et al. (2019) and apply a 1D convolutional derivative, which also avoids phase unwrapping.
In Appendix E, we show that the quad-polarimetric measurement (Fig. 2c) can be used to synthesize the full radar return from any antenna orientation using a matrix transformation

190
where γ is the angular offset from θ. Equation (18) is the mathematical equivalent to rotating the antennas in the field for each polarimetric configuration. As demonstrated in Fig. E1, we find no significant differences between the synthesized and the full azimuthal rotation dataset with 22.5°increments. Hence, if v3 is vertical a quad-polarimetric measurement is sufficient as it contains the full angular information.

Demonstration of anisotropic signatures in radar data using a synthetic model
For a given depth-profile of ∆λ(z), θ(z), and r(z), the radar return can be simulated using the forward model described by Eqs. (4)-(6). We show a seven layers synthetic model in Fig. 3 to visualize features in the radar data, which can be linked to ice-fabric parameters. The model parameters used to generate Fig. 3 are shown in Table 2. (yellow squares in blue areas), and v2 (yellow squares in red areas). The magnitude of Ψ at the black dots is the value of ∆λ.
Power anomalies illustrate the effects of anisotropic ice

200
where |s xx | is the amplitude of the complex received signal, and n is number of angular increments for θ. In P HH , a number of co-polarization nodes (CPN) occur, which result from destructive superposition of ordinary and extraordinary waves (Fig.   3a). The number of nodes per layer is only a function of ice fabric anisotropy in that layer, with higher horizontal anisotropy resulting in more nodes as predicted by Eq. (13). The nodes occur at a variable angular distance (termed AD in Fig. 3a) if anisotropic reflection is relevant (e.g., L2 and L3 in Fig. 3a). The angular dependency of the co-polarization nodes on 205 anisotropic scattering can be identified using a depth-invariant ice-fabric orientation (constant θ ). Previously, Fujita et al. (2006) approximated the correlation between AD and r with a linear regression. As detailed in Appendix F we improved this 12 https://doi.org/10.5194/tc-2020-370 Preprint. Discussion started: 20 January 2021 c Author(s) 2021. CC BY 4.0 License.
by finding the analytical solution Differences of both approaches are illustrated in Figure 4. Two important features in P HH are therefore the frequency of 210 occurrence of co-polarization nodes with depth (a first-order proxy for the horizontal anisotropy) and their angular distance (a mixed proxy for anisotropic reflections or depth-variable ice-fabric orientation). P HH can be 90°(e.g., L2) or 180°(e.g., L3) symmetric if r = 0 dB or r = 0 dB, respectively. This study Fujita et al. 2006 Figure 4. Dependence of reflection ratio on the azimuthal difference between two nodes as determined by Fujita et al. (2006) and through Eq. 20.
In a depth-invariant ice-fabric orientation, the minima in P HV align with v1 and v2 termed cross-polarization extinction (CPE in Fig. 3c). Using the radar forward model, this can be derived analytically for a single layer case as: In multi-layer cases, where θ changes with depth (e.g., L6 and L7 in Fig. 3b), P HV also depends on other parameters, making it difficult to infer θ using P HV alone.
The co-polarization nodes in P HH can also be observed in φ HHV V (termed DN in Fig. 3b). The depth of the node can 220 be automatically estimated at the zero-phase transition. Unlike P HH , the nodes in φ HHV V are 90°anti-symmetric, and their polarity is insensitive to r. This can be used to determine the directions of v1 and v2. The angular width of the nodes (termed AW in Fig. 3b) decreases when r = 0 dB (e.g., L3 or L4). The absolute value of Ψ at the principal axis's directions (v1 or v2) is a first-order proxy for ∆λ at a given depth (Eq. 17, Fig. 3d).
The vertical gridding of the model is 1 m. 3.5 An inverse approach to infer ice fabric from quad-polarimetric returns Our approach involves data preprocessing, initializing the model parameters, and parameter optimization using a constrained multivariable non-linear least-square inverse approach. All the three Eigenvalues are then estimated from the optimized ∆λ and r using a top to bottom layer-by-layer approach assuming isotropic ice at the surface.

Data preprocessing 235
The full angular response is synthesized from HH, VV, and HV observations for a single TR orientation (θ) using Eq. (18) at 1°increments. The amount and method of smoothing data depend on nodes' vertical frequency and phase polarity's sharpness.
The power anomalies are smoothed by moving average and 2D Gaussian convolution. The coherence phase (φ HHV V ) is inherently smoothed, depending on the size of the depth window in Eq. (14), while its gradient (Ψ) is smoothed with a 1D Gaussian convolution at each azimuth.

Model parameterization
We investigate two parameterization types for the free model parameters (θ, ∆λ, r) with depth: piece-wise constant and a superposition of Legendre Polynomials. The former has the highest number of free model parameters but can capture abrupt variability with depth. The latter has a reduced set of free model parameters with improved performance during the inversion, but varies more smoothly with depth. At Dome C, no abrupt variability is visible in the data so that we use the Legendre Polynomials with 50 free model parameters (30 for θ, 10 for ∆λ, and 10 for r). At EDML, we default to the piecewise constant parameterization resulting in 120 free model parameters (40 piecewise constant intervals at 50 m spacings for each model parameter).

Derivation of initial guess
The non-linear optimization problem depends on a well-defined initial guess based on our inferences from the synthetic data.

250
Initial guesses of variables are marked with superscript 0. We first derive the initial guess for the orientation of the ice fabric (θ 0 (z)) using the minima in P HV , polarity in φ HHV V , and the sign of Ψ. We then infer ∆λ 0 (z) using the absolute value of Ψ at the minima of P HV . The initial guess for r 0 (z) is zero. The underlying assumption for all of the initial guesses is that θ does not vary significantly with depth.

255
We optimize ∆λ, θ, and r for all depth intervals. There are a number of possible model data misfit metrics of power anomalies and phase differences 260 and the total misfit between the observed (obs.) and the modeled data (mod.) is defined as: where l 1 , l 2 , and l 3 are constants. In Table 3, we show the values of l 1 , l 2 , and l 3 that we used for Dome C and EDML sites. Using l 1 (for coherence phase misfit) in EDML is not applicable because the phase is too noisy. Therefore, we did not optimize ∆λ for this site. To further constrain the inversion, we set bounds on the model parameters so that 0 < ∆λ i < 0.5,

265
0°< θ i < 180°, and −30 dB < r i < 30 dB. This is implement in the cost function in the form of log-barrier functions using Matlab ® 's fmincon algorithm.

Reconstruction of all Eigenvalues
Once the radar forward model is optimized, we attempt to reconstruct all the three Eigenvalues in a top-to-bottom approach by solving Eq. 12 using the optimized r and ∆λ. At the surface ice is isotropic so that λ1 1 ≈ 0.33 allowing to infer λ2 1 and λ3 1 270 from ∆λ 1 For deeper layers i + 1, all three Eigenvalues, can in theory, be reconstructed analytically by solving

275
for λ1 i+1 and infer λ2 i+1 and λ3 i+1 with However, errors during the optimization may result in a reconstruction of the three Eigenvalues, which do not comply with limits inferred in Sect. 3.1. In that case, ∆λ and r are varied in a systematic search to find Eigenvalues within the permissible 280 limits. Solutions, in this case, are not unique, and additional constraints on the vertical gradients (here: 1.0 · 10 −6 < dλ3 dz < 1.5 · 10 −3 ) are required. This correction does not significantly alter the results from the previous section but assures that the inferred Eigenvalues are internally consistent.  (Figs. 5a, b). The existence of only one pair of nodes over 2000 m indicates comparatively small horizontal ice anisotropy (i.e., low ∆λ) similar to what has been observed at Dome Fuji (Fujita et al., 2006). The angular distance between the two co-polarization nodes is close to 90°, consistent with r close to 0 dB (Fig. 5a). P HV shows little depth-variability (Fig. 5c), suggesting that the ice-fabric orientation angle (θ) does not vary strongly with 290 depth. The scaled phase derivative (Ψ, Fig. 5d) is unclear in terms of polarity for the top 150 m. Below that, the polarity more clearly indicates the orientation of the largest horizontal Eigenvectors.
Optimized model results in Figs. 5e-h reproduce the principal patterns of the radar observations. The reconstructed Eigenvalues (Fig. 5i) capture the observed transition from isotropic to a girdle-type ice-fabric in the ice-core data. The reconstructed horizontal anisotropy (Fig. 5j)  observations. The ice-fabric orientation at the top 150 m is poorly constrained due to the low horizontal anisotropy (Fig. 5k).
The mean orientation of v2 below 150 m is 124°relative to True North, which is almost perpendicular to the surface flow direction towards 45°. The orientation cannot be validated with ice-core data, which is azimuthally unconstrained. The mean estimated reflection ratio below 150 m is low (r (z>150m) = −3 dB, Fig. 5l), indicating that the role of anisotropic reflections is small.

Ice-fabric parameters from polarimetric ApRES at EDML
Next, we apply our method to ApRES data collected at the EDML drill site. Contrary to what has been observed at Dome C, co-polarization nodes can barely be localized in P HH as no 90°symmetry is apparent (Fig. 6a). This indicates that anisotropic scattering is relevant (r = 0 dB), as already noticed earlier (Drews et al., 2012). Moreover, the coherence phase shows many nodes (Fig. 6b), indicating a much stronger horizontal anisotropy (i.e., large ∆λ). This is comparable to the ice core at Mizuho, 305 equally located in a flank flow regime (Fujita et al., 2006). Although P HV shows almost no depth variability in ice-fabric orientation (Fig. 6c), it is not straightforward to identify the direction of v1 and v2 using the polarity of Ψ because of the strong ice anisotropy (Fig. 6d).

The optimized model (Figs. 6e-h) reproduces all basic features seen in the radar data. Inferred model parameters closely
follow the ice-core measurements both in terms of absolute Eigenvalues (Fig. 6i) and horizontal anisotropy (Fig. 6j). The

Spatial variability of ice-fabric parameters in the local dome-flank transition zone
After investigating specific characteristics of a dome position (EDC) and a flank flow regime (EDML), we next investigate a local dome-to-flank transition (36 km). At Dome C, 19 sites are located along a profile extending 18 km away to either side 320 from the local ice dome (Fig. 1a), and a summary of the results is presented in Fig. 7. We focus on the upper 2000 m, where the signal to noise ratio is sufficiently high. All stations yield coherent results showing an isotropic ice fabric that gradually evolves into a weak girdle with depth. The depths of the first co-polarization nodes can be detected at all sites (green-dashed line in Fig. 7b). It is shallowest beneath the dome and moves to larger depths further away from the dome in the flanks. The depth-variability of the co-polarization nodes results in a ∆λ that is most anisotropic beneath the dome, and less anisotropic Our method extracts the horizontal anisotropy, orientation, and the anisotropic reflection ratio of the ice fabric as a function of 330 depth. We also estimate all three Eigenvalues required for the dielectric orientation tensor. This can be compared with laboratory measurements (Durand et al., 2009) and integrated into ice-flow models . Our main assumption is that the strongest Eigenvector (and with it the orientation tensor) is aligned in the vertical. We now discuss the limitations of our approach.
In terms of the data pre-processing, there are no structural differences in our data between synthesizing the polarization 335 dependency out of a single set of quad-polarimetric measurement (Appendix E) and the more common polarimetric measurements in glaciology where antennas are kept parallel or perpendicular while being rotated several increments between 0 and 180 degrees (Fujita et al., 2006). However, more systematic investigation is warranted if this also holds when the ice fabric is not aligned vertically.
The signal quality and noise level, particularly in the HHVV coherence phase, are important. In areas with high horizontal 340 anisotropy and consequently densely spaced co-and cross-polarization nodes (i.e., the EDML case), care needs to be taken that the denoising does not average over multiple nodes. Derivation of the initial guess for the inverse approach depends on the data quality and is guided by characteristic features in synthetic forward models, some of which can be analytically described for one layer cases. Multi-layer cases, however, are difficult to interpret, particularly if the ice-fabric orientation changes strongly (by several 10s of degrees) with depth. Fortunately, this does not appear to be the case for the data presented here, so that the 345 initial guess already results in a forward model that adequately captures characteristic features in the data. The optimization improves the model-data misfit but does not lead to significant differences with our first informed guess. Nevertheless, this step is required to predict the depth-variability of all the three Eigenvalues.
The reconstruction of the Eigenvalues assumes isotropic ice/firn at the surface. This is reasonable for the dome and flank-flow settings considered here, but may need to be revisited in other settings where fabric can develop near the surface as ice-streams 350 and ice-shelves. More critical is the reflection ratio itself, which is ill-constrained in magnitude and amplifies small changes in the Eigenvalues across the reflection boundaries. This is mitigated by the range of allowed Eigenvalues (Sect. 3.1), and it is those constraints that facilitate the derivation of all Eigenvalues from the anisotropic reflection ratio. The predicted Eigenvalues (λ1, λ2, and λ3) in this method such show a good match to the ice-core observations in both cases.
The azimuthal constraints that radar polarimetry provides can, in general, not be validated by ice-core measurements with 355 few exceptions (e.g., Westhoff et al., 2020). However, the alignment of the ice-fabric principal axes with the surface-flow direction detailed below adds credibility to our inferences and shows advantages of this approach over previous attempts focusing on the power anomalies only (Fujita et al., 2006;Matsuoka et al., 2012). The underlying reason for this is that the polarity of the depth gradient of the polarimetric coherence phase is independent of anisotropic scattering.
The inversion requires an initial guess (Sect. 3.5.3) that is based on experience from synthetic test cases. In our experience 360 with radar polarimetry and the explored ice dynamic context, this grants a robust solution, also because a wrong initial guess results in a large model-data misfit that can be identified easily. In the future, this can be improved by using gradient-free optimization schemes (e.g., in a Bayesian framework) that can correct for a poor initial guess by exploring the parameter space more systematically.
Our strongest assumption is that the strongest Eigenvector (v3) should be close to vertical. While this assumption is justified 365 here, as flow at domes is dominated by vertical compression and crystal c-axis tend to align in vertical, it may not apply elsewhere in ice sheets and cause an additional source of horizontal birefringence (Matsuoka et al., 2009). While it is possible to explore the effects of other than the largest Eigenvector being vertical (Jordan et al., 2019, p. 13), it is impossible to circumvent that the radio-wave propagation is vertical and hence insensitive to changes along that direction. In the future, we envision the use of wide-angle surveys with curved ray paths (e.g., Winebrenner et al., 2003) to overcome this limitation.

Spatial variability of ice-fabric types in a dome-flank transitions
We now turn to the ice-dynamic context of our inferred characteristics of ice-fabric variation from dome, where flow is dominated by vertical compression to flank flow regimes, where flow is driven by vertical shear. Our inverse approach shows higher horizontal ice anisotropy at EDML compared to Dome C throughout the ice column. This increase from the dome to the flank supports earlier inferences that ice anisotropy is larger in areas with significant horizontal strain compared to settings where 375 vertical compression is dominant (Fujita et al., 2006;Matsuoka et al., 2012). This is in contrast, however, with the observed decrease in ice anisotropy in the Dome C transect (Fig. 7c), where the ice fabric is more anisotropic at the Dome compared to the flanks. Our hypothesis is that this near-field anomaly reflects ice-dynamic modification of ice fabric through the Raymond effect (Raymond, 1983). Martín et al. (2009)  hence approximately covers this domain. The absence of Raymond arches in the radar stratigraphy beneath Dome-C (Cavitte et al., 2016, p. 325) suggests that these need a longer time to evolve, whereas the ice-fabric pattern reflects the instantaneous operation of the Raymond effect. We acknowledge that there are other explanations for the ice-fabric pattern under Dome C, such as across-profile flow or bedrock influence. In any case, we want to highlight here how, due to the spatial extension of our observations, our inferred fabric distributions combined with an anisotropic flow model can be used to test these and other 385 hypothesis.
Both in EDML and Dome C areas, the inferred ice-fabric orientation varies little over the depth-intervals considered, and in both cases, the inferred orientations line-up with the surface flowfield. More specifically, v1 is approximately oriented along-flow and v2 is approximately oriented across-flow. Those directions also align with the principal strain rate components in Dome C (Rémy and Tabacco, 2000;Vittuari et al., 2004) and EDML (Drews et al., 2012) (Fig. 1). In both cases, v2 is factor to explain the radar signatures, while at Dome C, it is close to negligible. Fujita et al. (2006) have observed a similar increase in anisotropic scattering between Dome-F and Mizuho, suggesting that this may be a generic feature in ice sheets that requires more investigation. Contrary to EDML, the signal at Dome C is dominated by birefringence, and the contribution of anisotropic reflection is small. Yet, it appears that it leaves a small signature in the data that can be exploited. Moreover, our analysis suggests that there are no other mechanisms (e.g., a directional interface roughness) contributing to anisotropic 400 reflections. This point requires confirmation from other ice-core sites because the recovery of all three Eigenvalues (and their corresponding directions) offers significant possibilities to constrain ice fabric in ice sheets in general.

Conclusions
We show here, for the first time, the spatial distribution of ice-fabric in domes: from the summit, where flow is dominated by vertical compression, to the flanks, where flow is driven by vertical shear. The combination of co-and cross-polarized power 405 anomaly along with the depth gradient of polarimetric coherence phase provides three major parameters and their changes over depth, i.e., the ice-fabric orientation, horizontal anisotropy, and its vertical variability. We quantify these changes using an inverse approach that extracts ice-fabric information from radar polarimetry. We present here a method to combine them and infer the full orientation tensor. We validate our technique with data from two ice-core locations situated in contrasting ice-flow regimes. The inferred ice-fabric orientation aligns with the observed surface velocity and surface strain rate fields.

410
This suggests that polarimetric radar is an ideal tool to map ice-fabric characteristics elsewhere as well.
We present ice-fabric spatial distribution across a flow-plane at Dome C. The 20 ApRES sites in that area are internally consistent, and small changes in the horizontal anisotropy can horizontally be tracked in the polarimetric coherence phase. We detect a minor decrease in horizontal anisotropy away from the dome that we tentatively link to the operation of the Raymond effect. On larger spatial scales, the horizontal anisotropy increases in the flanks (i.e., at EDML), and our findings are consistent 415 with previous studies. Our analysis suggests that ice-fabric characteristics can now be reliably inferred in larger parts of the Antarctica and Greenland ice sheet, given that more and more profiles are recorded in coherent and in the quad-polarimetric configuration. This will be a decisive step to further constrain the anisotropic nature of ice and understand better its contribution to internal deformation.

Appendix D: Polarity of the coherence phase gradient
This section details the relationship between the polarity of the phase gradient and the corresponding directions of the Eigenvectors. Care has to be taken here, as the de-ramping during ApRES data acquisition is equivalent to a complex conjugation of the received signal. If this is not accounted for, the inferred Eigenvector v1 and v2 will be swapped. More specifically, for a 450 received signal at θ = 0°: s V V = A (Γ y cos(2zk y ) + jΓ y sin(2zk y )) , so that the coherence phase results in: and conversely for θ = 90°: As explained in Sect. 3.3, k x and k y are a function of λ1 and λ2, respectively. Because λ1 ≤ λ2 it follows that k x < k y .
Therefore, φ HHV V (θ = 0°) < 0 and φ HHV V (θ=0°) dz < 0. The reverse holds for θ = 90°. Appendix E: Reconstruction of azimuthal measurements from a single quad-polarimetric acquisition The transformation is purely geometrical and corresponds to a coordinate transformation into a rotated reference system for an arbitrary γ: resulting in: The remaining quadratic equation has two solutions corresponding to the two co-polarization nodes: The angular distance between these nodes then results in which can be re-arranged for the reflection ratio as: r = 1 tan 2 ( AD 2 ) .
contributed to the writing of the final manuscript.
Competing interests. OE is Co-Editor in Chief and RD is Editor of The Cryosphere.
Acknowledgements. Acknowledgements: RE and RD were supported by a DFG Emmy Noether grant DR 822/3-1. This publication was also generated in the frame of Beyond EPICA. The project has received funding from the European Union's Horizon 2020 research and innovation programme under grant agreement No. 815384 (Oldest Ice Core) and 730258 (CSA). It is supported by national partners and 490 funding agencies in Belgium, Denmark, France, Germany, Italy, Norway, Sweden, Switzerland, The Netherlands and the United Kingdom.
The Dome C measurements were made possible by the logistic provided by IPEV (prog. 902) and PNRA. We thank Luca Vittuari (University of Bologna, Italy) for the positioning of the stakes. The opinions expressed and arguments employed herein do not necessarily reflect the official views of the European Union funding agency or other national funding bodies. This is Beyond EPICA publication number XX. We thank the AWI logistics personnel for support of the work at Kohnen.