Gravity Topography Modeling of the Denman Glacier Region Using a Geostatistical Approach

Lösing, Mareen; Aitken, Alan; Field, Michael; MacKie, Emma; Li, Lu

doi:10.5194/egusphere-2026-166

Preprints

https://doi.org/10.5194/egusphere-2026-166

Preprints

10 Feb 2026

| 10 Feb 2026

Gravity Topography Modeling of the Denman Glacier Region Using a Geostatistical Approach

Mareen Lösing, Alan Aitken, Michael Field, Emma MacKie, and Lu Li

Abstract. The Denman Glacier is one of East Antarctica’s most dynamic outlet systems and is modeled to host the deepest continental marine trough, with the potential to contribute up to ∼ 1.5 m of global sea level rise. Yet its bed geometry remains poorly constrained because airborne radar surveys struggle to image steep, narrow, and deeply incised troughs, and existing continental-scale compilations rely heavily on interpolation or mass-conservation assumptions. During the Australian Denman Terrestrial Campaign 2023/24, high-resolution ground-based gravity measurements were collected across the deepest part of the trough, providing short-wavelength constraints that complement ICECAP airborne gravity and radar data.

We use a two-scale, ensemble-based gravity inversion to reconstruct Denman’s bed topography. A geostatistical separation of terrain and non-terrain gravity effects generates multiple plausible regional background fields, and each is explored with a random-walk Metropolis–Hastings Markov Chain Monte Carlo (MCMC) inversion in which small, spatially correlated Gaussian perturbations modify the bed. Within each realization, candidate geometries are jointly evaluated against the gravity signal and radar picks, with gravity providing the dominant constraint in regions lacking radar coverage.

The resulting ensemble reveals a more rugged, spatially variable, and internally segmented subglacial landscape than represented in current bed products. Along cross-trough profiles, the best-fit gravity-derived bed generally falls between BedMachine and Bedmap3 yet exhibits steeper trough walls and greater lateral relief. Depth-to-magnetic-source estimates further support contrasting lithologies across the trough, with crystalline basement to the west and weaker, sedimentary signatures to the east.

The inferred geometry, characterized by steep flanks and an inland-sloping basin, reinforces the susceptibility of Denman Glacier to Marine Ice Sheet Instability. These results highlight the importance of incorporating geophysical inversion methods into future Antarctic bed-mapping efforts and provide an ensemble of bed realizations suitable for ice flow modeling and assessments of grounding line stability.

Received: 13 Jan 2026 – Discussion started: 10 Feb 2026

Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors. Views expressed in the text are those of the authors and do not necessarily reflect the views of the publisher.

Download & links

Mareen Lösing, Alan Aitken, Michael Field, Emma MacKie, and Lu Li

Status: final response (author comments only)

RC1:
'Comment on egusphere-2026-166', Matthew Tankersley, 27 Mar 2026
Lösing et al. used a novel gravity inversion scheme to model to bed topography beneath Antarctica’s Denman Glacier, revealing a steep-flanked and inland-sloping narrow trough. Their approach utilized existing airborne gravity data, and a recently collected ground-based gravity profile. The geostatistical modelling approach, compared to traditional deterministic gravity inversions, provides an avenue for the much-needed estimates of spatial variability in model uncertainties. While I think there is some work to be done to clarify some of the methods and approaches, I think this paper makes a substantial and important contribution to the field.

General comments:
- For the isostatic anomaly calculation, you cite Liebsch 2020, but I could only find some code associated with this. Can you add 1 or 2 sentences in the text or appendix describing this method? Including which bed topography model you used, and which density values.

- On line 174 you state you will use gravity disturbances instead of FAA. I would change all further references of FAA to gravity disturbance, such as in Fig 2 and its caption, and the equations in sections 3.1. Also I would initially state that the terrain effects are the same thing as the Bouguer correction, and that non-terrain effects are the same thing as the Bouguer disturbance, but once you state this, stick with one naming convention.

- I’ve added lots of suggestions for the text of Steps 1-6. I think this is part is especially important for readers to follow, and I struggled to understand exactly what you did. Hopefully my suggestions can help clarify the workflow you used and help compare/contrast it to the methods of other bathymetry inversion studies to put your methodological advance in context.

- Clarification of gravity anomalies:
From the code you provided, this seems to be your workflow:
for clarity, bouguer correction = terrain effects, and bouguer disturbance == non-terrain effects
1) subtract isostatic gravity effect from FAA to get an isostatically-corrected FAA
2) calculate the bouguer correction of ice, water, and rock from BedMachine and subtract this from the isostatically-corrected FAA to get the bouguer disturbance
3) calculate a smoothed (regional) bouguer disturbance with a RBF of the data only near constraints
4) Subtract the regional bouguer disturbance from the full bouguer disturbance to get a residual bouguer disturbance.
- up to this point, this should be equivalent to An et al. 2015 DC shift or Tankersley et al. 2025 constraint point regionalization
5) cut out the residual bouguer disturbance away from radar data, and replace it with SGS to get an ensemble of residual bouguer disturbances.
6) add the regional bouguer disturbance back to each of the residual bouguer disturbance ensemble to get an ensemble of full bouguer disturbances
7) subtract this ensemble of full bouguer disturbances from the isostatically-corrected FAA to get an ensemble of bouguer corrections

I think this differs slightly from the presentation in Figure 2, which omits the important step of calculating the regional bouguer disturbance, and makes it seem like the SGS is applied directly to the full bouguer disturbance, instead of to the residual component.

I my interpretation of the code / your methodology is correct, I would suggest somewhere you show all the following maps of gravity data, either in figure 3 or a supplementary figure: FAA (Fig 3a right?), moho effect, isostatically-corrected FAA, the mean terrain effect, the mean Bouguer disturbance, and the regional and residual components of the Bouguer disturbance. I think Fig 3. would be more informative if it showed 1) the isostatically-corrected FAA (instead of the FAA), 2) the mean target terrain effect (current Fig 3b), and the mean gravity misfit (forward effect of Bedmachine - mean target terrain effect). I would then make sure to cross-reference which subplots of Fig 2 correspond to which subplots of Fig 3 / the new supplementary figure.

- Corrections for ice and water: I read through your code trying to understand how you correct for the gravity effect of ice and water, specifically your python class "PrismGen2". It seems you make 6 sets of prisms, a positive and negative set each for ice, water, and rock, but it was pretty hard to follow. Could you clarify a bit in the text, or supplement, about how you discretized and what density values you used for each of these 6 prism sets, and whether all 6 prism sets were computed at each iteration (or just the bed prisms?).

Specific comments:
Line 44: You might add (Fig. 1c) after "These gaps".
Line 98: Maybe use ground survey instead of ground line to avoid confusion with grounding line.
Line 117: "above the ice surface" suggests these altitudes are with respect to the ice surface. Are they? Or is it 1.5-5 km ellipsoidal heights?
Line 131: Did you block-reduce the line data before fitting the equivalent sources? If so, what cell size?
Line 137: Was the block median here used to go from the 2km resolution gravity grid to the 1km? Why not just predict the equivalent sources at both the 1 and 2km grids? Or does block-reduction essentially give the same result?
Line 138: The fact that you mask the gravity data to 2km is quite novel, as most inversions use full grids without any holes. If you retain the inverted bed solution beyond 2km from gravity data, why limit the gravity to this distance? How would you expect the results to change if you used a fully-interpolated grid of gravity data?
Line 181: Is the ground-gravity data already at the surface elevation, since that is where they were collected? So the AntGG data was only upward continued for the ICECAP comparison, since the AntGG data was already at the surface? Or was the ground-gravity data upward continued to 3.5 km as well?
Line 188: "geometry" is a little confusing, maybe "bed topography"?
Line 188: since it's the first mention of it, maybe say "non-terrain gravity effects" for clarity.
Line 188: add "and remove" after isolate
Line 197: I think this could be re-worded to be clearer, for example: We therefore (i) forward-model and subtract the ‘known’ terrain effect from the isostatically-corrected gravity disturbance, giving a Bouguer disturbance, (ii) estimate a regional component of the Bouguer disturbance by fitting a smooth trend to the gravity data at carefully selected observations points with topography constraints, (iii) remove this regional component to obtain the residual Bouguer disturbance, and (iv) at locations far from topography constraints, replace these residual values with an ensemble of simulate values.
Line 200: It might be helpful to add a reference to Fig. 2A for this section, as that is what it is describing.
Line 213: By “remaining broad signals not explained by mapped topography”, are you referring to the gravity effect of the bed resulting from differences between the true bed and Bedmachine? If so, I would replace the word broad with residual, as these signals are not broad, as there are only short wavelength deviations between the true topography and bedmachine topography. I think adding the below sentence will help clarify: This residual Bouguer disturbance is our desired input to the inversion, so we first need to estimate and remove the regional component.
Line 216: I suggest changing the sentence to as follows for clarity. "It consists of all bouguer disturbance grid points outside the inversion domain and, within the domain, only bouguer disturbance grid points co-located with radar picks ...."
Line 218: This section is a little confusing, as it's discussing a conditional gravity subset, but ends with saying "Interior points ... treated as unknown topography", which seems to refer to a topography data subset, not a gravity subset. Maybe just remove “and treated as unknown topography.”
Line 220: As above, I think a little clarification could help the reader. "broad-scale trend" -> "the regional (or can say long-wavelength) component of the gravity disturbance"
Line 225: Up to this point, this technique is the same as the DC shift technique used in An et al. 2015, or the constraint-point-minimization used in Tankersley et al 2025, just with a difference of the interpolation algorithm (radial-basis function, vs min.curvature / biharmonic splines), correct?. This might be good to explicitly state to help put your method in context of other literature.
Line 232: short-wavelength residual -> short-wavelength Bouguer disturbance residual
Line 233: Add reference to Fig. 2b center.
Line 238: Add reference to Fig. 2b right side.
Line 240: This states B is the sum of T (the regional bouguer disturbance) and r (the residual bouguer disturbance), so I would remove the word regional before “non-terrain effects”, since it’s not the regional component, but the entire bouguer disturbance.
Line 243: observations -> isostatically-corrected gravity disturbance
Line 308: the gridcell the radar data falls in is replaced with the radar data elevation. Does this not lead to pedestals or holes? How do you achieve a smooth transition? Do you do this for the actually bed inversion depths, or just for computing the density model?
Line 323: the DOI link for Melo and Barbosa is incorrect
Line 336: I think you mean to refer to fig 5c.
Line 414: grounding line zone -> grounding zone, same on line 424 and caption of Fig 10.

Fig 2.
"Cutout gravity away from known topography" might be clearer than "Cutout unknown Topography", since you are cutting out the gravity effect, not the topography, right?

I think this figure should more clearly relate to your steps 1-6 in the text. For example, it is unclear where, or if, step 3 (regional separation) is shown.

Fig 3.
a) Is the FAA in ‘a’ adjusted for the Moho gravity effect already? If so, I would change the title to “Isostatically-corrected gravity disturbance”

b) change the title to target terrain effect, and optionally leave Mean g_target in parenthesis

c) change title to bouguer disturbance (can leave mean B(x) in parentheses)

Fig 5.
Can you replace “major geographic boundaries” with the grounding line, or is it different? Same in Fig. 6 and 7.

Fig.7. Could you add a small histogram of the residuals for each of these figures? Then can be stacked vertically to the right side of the figures.

Fig 10. The black line is showing the grounding line not the coast line right?

Fig A1. I think adding to the caption that the right side of the x-axis corresponds to ~ 58 days is useful.
Citation: https://doi.org/10.5194/egusphere-2026-166-RC1
- AC1: 'Reply on RC1', Mareen Lösing, 15 Jun 2026
  
  Thank you very much for your helpful comments. We agree that the terminology used in the original manuscript was not always consistent and may have caused confusion. Following your suggestion, we have revised the text to consistently distinguish between gravity disturbances, isostatically corrected gravity disturbances, terrain effects, and non-terrain gravity disturbances throughout the manuscript.
  To improve clarity, we have also expanded the description of the gravity-field decomposition workflow and added a dedicated subsection describing the isostatic correction. This section now explains the Airy-compensated Moho model, the density assumptions, and the calculation of the resulting isostatically corrected gravity disturbance. We further revised Figures 2 and 3, their captions, and the corresponding equations to ensure that the terminology is used consistently and that the relationship between the different gravity-field components is explicit.
  We believe these changes substantially improve the readability of the manuscript and reduce potential ambiguity regarding the gravity quantities used in the inversion workflow.
  
  Line-by-line responses:
  - For the isostatic anomaly calculation, you cite Liebsch 2020, but I could only find some code associated with this. Can you add 1 or 2 sentences in the text or appendix describing this method? Including which bed topography model you used, and which density values.
  We added a separate subsection (Section 2.4) to better explain the isostatic correction, including the underlying method, topographic model, and density assumptions.
  - On line 174 you state you will use gravity disturbances instead of FAA. I would change all further references of FAA to gravity disturbance, such as in Fig 2 and its caption, and the equations in sections 3.1. Also I would initially state that the terrain effects are the same thing as the Bouguer correction, and that non-terrain effects are the same thing as the Bouguer disturbance, but once you state this, stick with one naming convention.
  We added a clarification on the naming of terrain effects and non-terrain disturbances and removed "FAA" from figures, captions, and equations.
  - I’ve added lots of suggestions for the text of Steps 1-6. I think this is part is especially important for readers to follow, and I struggled to understand exactly what you did. Hopefully my suggestions can help clarify the workflow you used and help compare/contrast it to the methods of other bathymetry inversion studies to put your methodological advance in context.
  
  - Clarification of gravity anomalies:
  
  From the code you provided, this seems to be your workflow:
  
  for clarity, bouguer correction = terrain effects, and bouguer disturbance == non-terrain effects
  
  1) subtract isostatic gravity effect from FAA to get an isostatically-corrected FAA
  
  2) calculate the bouguer correction of ice, water, and rock from BedMachine and subtract this from the isostatically-corrected FAA to get the bouguer disturbance
  
  3) calculate a smoothed (regional) bouguer disturbance with a RBF of the data only near constraints
  
  4) Subtract the regional bouguer disturbance from the full bouguer disturbance to get a residual bouguer disturbance.
  
  - up to this point, this should be equivalent to An et al. 2015 DC shift or Tankersley et al. 2025 constraint point regionalization
  We added these references in the Methods section to make this conceptual similarity clearer.
  5) cut out the residual bouguer disturbance away from radar data, and replace it with SGS to get an ensemble of residual bouguer disturbances.
  
  6) add the regional bouguer disturbance back to each of the residual bouguer disturbance ensemble to get an ensemble of full bouguer disturbances
  
  7) subtract this ensemble of full bouguer disturbances from the isostatically-corrected FAA to get an ensemble of bouguer corrections
  
  I think this differs slightly from the presentation in Figure 2, which omits the important step of calculating the regional bouguer disturbance, and makes it seem like the SGS is applied directly to the full bouguer disturbance, instead of to the residual component.
  
  I my interpretation of the code / your methodology is correct, I would suggest somewhere you show all the following maps of gravity data, either in figure 3 or a supplementary figure: FAA (Fig 3a right?), moho effect, isostatically-corrected FAA, the mean terrain effect, the mean Bouguer disturbance, and the regional and residual components of the Bouguer disturbance. I think Fig 3. would be more informative if it showed 1) the isostatically-corrected FAA (instead of the FAA), 2) the mean target terrain effect (current Fig 3b), and the mean gravity misfit (forward effect of Bedmachine - mean target terrain effect). I would then make sure to cross-reference which subplots of Fig 2 correspond to which subplots of Fig 3 / the new supplementary figure.
  Thank you for this suggestion. We have substantially revised both Figures 2 and 3 following your recommendation. We have also updated Section 3.1 and added explicit cross-references in the text to the relevant panels in Figures 2 and 3 to clarify the relationship between the gravity-field decomposition and the inversion inputs. We believe these changes make the workflow and interpretation considerably clearer.
  - Corrections for ice and water: I read through your code trying to understand how you correct for the gravity effect of ice and water, specifically your python class "PrismGen2". It seems you make 6 sets of prisms, a positive and negative set each for ice, water, and rock, but it was pretty hard to follow. Could you clarify a bit in the text, or supplement, about how you discretized and what density values you used for each of these 6 prism sets, and whether all 6 prism sets were computed at each iteration (or just the bed prisms?).
  In the updated Section, we now explicitly describe how the model domain is discretized and clarify the densities used.
  Specific comments:
  Line 44: You might add (Fig. 1c) after "These gaps".
  Added.
  Line 98: Maybe use ground survey instead of ground line to avoid confusion with grounding line.
  Changed.
  Line 117: "above the ice surface" suggests these altitudes are with respect to the ice surface. Are they? Or is it 1.5-5 km ellipsoidal heights?
  
  We have clarified this in the manuscript as follows:
  
  “Surveys were conducted at varying aircraft heights, typically ranging from approximately 1.5 km to ~4 km ellipsoidal heights.”
  Line 131: Did you block-reduce the line data before fitting the equivalent sources? If so, what cell size?
  The equivalent sources were fitted directly to the original line data without prior block reduction. The upward-continued field was then predicted onto the BedMachine reference grid in 1 and 2 km resolution for local and regional inversion, respectively.
  Line 137: Was the block median here used to go from the 2km resolution gravity grid to the 1km? Why not just predict the equivalent sources at both the 1 and 2km grids? Or does block-reduction essentially give the same result?
  Thank you for pointing this out. The wording was unclear. The block-median reduction was not applied to the upward-continued gravity field. Instead, it was used to resample the BedMachine topography onto the 1 km and 2 km reference grids used for the local and regional inversions. The equivalent-source model was then evaluated directly at the nodes of these predefined grids. We have revised the text to clarify this workflow.
  Line 138: The fact that you mask the gravity data to 2km is quite novel, as most inversions use full grids without any holes. If you retain the inverted bed solution beyond 2km from gravity data, why limit the gravity to this distance? How would you expect the results to change if you used a fully-interpolated grid of gravity data?
  We have revised the manuscript to better clarify the rationale for the 2 km masking. The mask is chosen to ensure that only gravity values directly supported by observations are used, avoiding artifacts from equivalent-source extrapolation. This reflects a deliberately conservative approach, where we want to avoid forcing the inversion to fit values in regions without observational constraint. Using a fully interpolated gravity grid would likely produce a smoother and more spatially continuous result.
  Line 181: Is the ground-gravity data already at the surface elevation, since that is where they were collected? So the AntGG data was only upward continued for the ICECAP comparison, since the AntGG data was already at the surface? Or was the ground-gravity data upward continued to 3.5 km as well?
  The ground-based gravity data were collected at (near-)surface elevation and were not upward continued. Upward continuation was applied only to the AntGG2021 field for comparison with the ICECAP data. We have clarified this in the revised manuscript.
  Line 188: "geometry" is a little confusing, maybe "bed topography"?
  Changed.
  Line 188: since it's the first mention of it, maybe say "non-terrain gravity effects" for clarity.
  Changed.
  Line 188: add "and remove" after isolate
  Changed.
  Line 197: I think this could be re-worded to be clearer, for example: We therefore (i) forward-model and subtract the ‘known’ terrain effect from the isostatically-corrected gravity disturbance, giving a Bouguer disturbance, (ii) estimate a regional component of the Bouguer disturbance by fitting a smooth trend to the gravity data at carefully selected observations points with topography constraints, (iii) remove this regional component to obtain the residual Bouguer disturbance, and (iv) at locations far from topography constraints, replace these residual values with an ensemble of simulate values.
  Thanks for this comprehensive rewording. We have significantly updated this paragraph, and hope this will improve clarity of the workflow.
  Line 200: It might be helpful to add a reference to Fig. 2A for this section, as that is what it is describing.
  Added.
  Line 213: By “remaining broad signals not explained by mapped topography”, are you referring to the gravity effect of the bed resulting from differences between the true bed and Bedmachine? If so, I would replace the word broad with residual, as these signals are not broad, as there are only short wavelength deviations between the true topography and bedmachine topography. I think adding the below sentence will help clarify: This residual Bouguer disturbance is our desired input to the inversion, so we first need to estimate and remove the regional component.
  We removed the word ‘broad’.
  Line 216: I suggest changing the sentence to as follows for clarity. "It consists of all bouguer disturbance grid points outside the inversion domain and, within the domain, only bouguer disturbance grid points co-located with radar picks ...."
  Changed.
  Line 218: This section is a little confusing, as it's discussing a conditional gravity subset, but ends with saying "Interior points ... treated as unknown topography", which seems to refer to a topography data subset, not a gravity subset. Maybe just remove “and treated as unknown topography.”
  Changed.
  Line 220: As above, I think a little clarification could help the reader. "broad-scale trend" -> "the regional (or can say long-wavelength) component of the gravity disturbance"
  Changed.
  Line 225: Up to this point, this technique is the same as the DC shift technique used in An et al. 2015, or the constraint-point-minimization used in Tankersley et al 2025, just with a difference of the interpolation algorithm (radial-basis function, vs min.curvature / biharmonic splines), correct?. This might be good to explicitly state to help put your method in context of other literature.
  As mentioned above, we added these references to the Methods section to better place our approach in the context of previous gravity inversion studies.
  Line 232: short-wavelength residual -> short-wavelength Bouguer disturbance residual
  Changed.
  Line 233: Add reference to Fig. 2b center.
  Added.
  Line 238: Add reference to Fig. 2b right side.
  Added.
  Line 240: This states B is the sum of T (the regional bouguer disturbance) and r (the residual bouguer disturbance), so I would remove the word regional before “non-terrain effects”, since it’s not the regional component, but the entire bouguer disturbance.
  Changed.
  Line 243: observations -> isostatically-corrected gravity disturbance
  Changed.
  Line 308: the gridcell the radar data falls in is replaced with the radar data elevation. Does this not lead to pedestals or holes? How do you achieve a smooth transition? Do you do this for the actually bed inversion depths, or just for computing the density model?
  This replacement is applied only within the forward modelling step used to refine the density field, rather than modifying the inversion result itself. Sharp discontinuities are mitigated by the regularisation of the inversion and the spatial sensitivity of the gravity response. We updated this paragraph accordingly.
  Line 323: the DOI link for Melo and Barbosa is incorrect
  Changed.
  Line 336: I think you mean to refer to fig 5c.
  Changed.
  Line 414: grounding line zone -> grounding zone, same on line 424 and caption of Fig 10.
  Changed.
  Fig 2. "Cutout gravity away from known topography" might be clearer than "Cutout unknown Topography", since you are cutting out the gravity effect, not the topography, right?
  
  I think this figure should more clearly relate to your steps 1-6 in the text. For example, it is unclear where, or if, step 3 (regional separation) is shown.
  We have significantly revised Figure 2 to improve the clarity of the workflow and its correspondence with the methodological steps described in the text.
  Fig 3. a) Is the FAA in ‘a’ adjusted for the Moho gravity effect already? If so, I would change the title to “Isostatically-corrected gravity disturbance”
  
  b) change the title to target terrain effect, and optionally leave Mean g_target in parenthesis
  
  c) change title to bouguer disturbance (can leave mean B(x) in parentheses)
  We have changed this Figure significantly.
  Fig 5. Can you replace “major geographic boundaries” with the grounding line, or is it different? Same in Fig. 6 and 7.
  Thank you for this suggestion. The boundaries are derived from the BedMachine v3 mask and represent the extent of floating ice (mask = 3) and ice-free land (mask = 1). These contours therefore include the ice-shelf front and the land boundary, which corresponds to the line where land is exposed, but also to ice/land–ice-shelf interfaces.
  As a result, these boundaries do not strictly correspond to the grounding line. We have therefore retained the more general wording but clarified this explicitly in the figure captions to avoid confusion.
  Fig.7. Could you add a small histogram of the residuals for each of these figures? Then can be stacked vertically to the right side of the figures.
  We added histograms to Figures 7.
  Fig 10. The black line is showing the grounding line not the coast line right?
  As explained above, we have clarified this in the captions.
  Fig A1. I think adding to the caption that the right side of the x-axis corresponds to ~ 58 days is useful.
  Changed.
  
  Citation: https://doi.org/10.5194/egusphere-2026-166-AC1
RC2:
'Comment on egusphere-2026-166', Anonymous Referee #2, 06 May 2026

This paper sets out to invert for subglacial topography (and bathymetry) of the Denman Glacier region using airborne gravity data, augmented with local ground-based gravity observations. The method uses a geostatistical approach to create an ensemble of realistic terrain gravity effects theoretically isolated from the background geology. These terrain gravity effects are then inverted for topography using a Markov Chain Monte Carlo (MCMC) approach, retaining uncertainty associated with the regions where gravity data alone is inverted.
Overall, the paper is well written and the methods are clear. However, I have two key concerns/clarifications which I think need addressing before publication, as well as a few other minor comments.
First – how is the ground gravity data incorporated into the raster of gravity data used in the inversion? L132 states “we upward continued the data to 3.5 km”, and this data was interpolated onto 1 km and 2 km rasters for inversion. However, this section seems to just be referring to the airborne data. Fig 3a shows distinct negative values (if you really zoom in) which suggest the ground data was incorporated, but how this was done should be explained. It may be better to describe the ground gravity data 1^st then the airborne, including how the land data was added. Also – the ground and airborne data directly overlie one another. If a grid was created by a simple block-median reduction the oversampled airborne gravity data would significantly reduce the utility of the nine ground based measurements which do however theoretically have significantly higher resolution due to their lower altitude and without along track filtering.
Second – the method appears to generate a certain amount of noise with wavelengths of 5 to 10 km away from any observational data (gravity or radar). It is not clear what the origin of this noise is and the amplitude is hard to determine. However, it appears to be of a similar form to the Gaussian random perturbation shown in Fig. 4 b and d. My concern is that the method may be introducing noise (bumps and depressions) where there is no real data. It then becomes harder to be confident of the bumps in the bed of the main channel are real or inversion artifacts. One way to convince readers the recovered bed is really due to signal in the gravity data is to plot the gravity data above the topography (like in Fig. 9). Alternatively or in addition it would be good to quantify the amplitude of the noise imposed away from areas with direct constraint, as this could be factored in for example when considering the reliability of water routing across a region where teh bed was derived using this technique.
Line by line points:
L132 Is it only the airborne data which is upward continued, or the ground based data as well?
L175 – it is unclear how ground and airborne data are combined.
L184-186 It is not clear why a different correction was applied along the airborne profile corresponding to the ground data compared to the wider airborne dataset.
L205 What observation altitude was used in the terrain model? I presume it is the same as the previously noted upward continuation altitude, but this should be stated.
L294 I understand you repeat the MCMC approach for each ensemble member to estimate the error, but how do you choose the ‘best’ model, or is it taken as the ensemble average of all the recovered topographies? I think this is reported on L335, but not sure.
L372-374 Suggests the gravity inversion shows topography “characterised by localised depressions and subdued ridges, which are muted or absent in current continent-scale compilations”. This is true, but are these features seen in the gravity data or artifacts of the inversion? As noted above these undulations are seen away from any data (gravity or radar). It may be that the inversion is generating realistic looking texture. This could be good as the undulations provide realistic texture, but it must be acknowledged that away from gravity and radar observations this is still a form of artificial noise.
L376 MAE – this is the first occurrence of this abbreviation so should be defined.
L453-456 These lines are duplicated immediately below (L457-L460).

Citation: https://doi.org/10.5194/egusphere-2026-166-RC2
- AC2: 'Reply on RC2', Mareen Lösing, 15 Jun 2026
  
  Thank you for your positive assessment of the manuscript and for the thoughtful and constructive comments. We agree that small-scale bed features in regions lacking direct gravity or radar constraints should not be interpreted as uniquely resolved topography. In these areas, the short-wavelength component is generated through conditional stochastic simulations and therefore represents plausible morphological variability rather than deterministic bed structure.
  Line-by-line responses:
  First – how is the ground gravity data incorporated into the raster of gravity data used in the inversion? L132 states “we upward continued the data to 3.5 km”, and this data was interpolated onto 1 km and 2 km rasters for inversion. However, this section seems to just be referring to the airborne data. Fig 3a shows distinct negative values (if you really zoom in) which suggest the ground data was incorporated, but how this was done should be explained. It may be better to describe the ground gravity data 1st then the airborne, including how the land data was added. Also – the ground and airborne data directly overlie one another. If a grid was created by a simple block-median reduction the oversampled airborne gravity data would significantly reduce the utility of the nine ground based measurements which do however theoretically have significantly higher resolution due to their lower altitude and without along track filtering.
  Thank you for pointing this out. We agree that the previous description was unclear. We have added a paragraph in Section 2.2 to clarify that we added these data at their original heights explicitly after upward continuation of the airborne data.
  Second – the method appears to generate a certain amount of noise with wavelengths of 5 to 10 km away from any observational data (gravity or radar). It is not clear what the origin of this noise is and the amplitude is hard to determine. However, it appears to be of a similar form to the Gaussian random perturbation shown in Fig. 4 b and d. My concern is that the method may be introducing noise (bumps and depressions) where there is no real data. It then becomes harder to be confident of the bumps in the bed of the main channel are real or inversion artifacts. One way to convince readers the recovered bed is really due to signal in the gravity data is to plot the gravity data above the topography (like in Fig. 9). Alternatively or in addition it would be good to quantify the amplitude of the noise imposed away from areas with direct constraint, as this could be factored in for example when considering the reliability of water routing across a region where teh bed was derived using this technique.
  Thank you for this comment. We agree that in areas far from constraints, the short-wavelength bed features are not directly constrained by gravity. Instead, it is generated by Gaussian random fields. We have therefore clarified that these areas should be interpreted cautiously. At the same time, the standard deviation of the model ensemble (Fig 5a) shows a general trend with an increase in model uncertainty further away from data constraints. However, areas with gravity constraints without radar constraints have the highest uncertainties (e.g. under the main trough and under ice-shelf). It suggests the nature of non-uniqueness of the gravity inversion, and also the interpolated method might underestimate the bed uncertainties. In areas without gravity constraints, our method and the smooth interpolation method are all conditioned by radar measurements. The ensemble-based approach allows for more geologically realistic spatial variability and morphological information than would be obtained from purely smooth interpolation products.
  Line by line points:
  L132 Is it only the airborne data which is upward continued, or the ground based data as well?
  It is only the airborne data. We have clarified this more.
  L175 – it is unclear how ground and airborne data are combined.
  As in the comment above, we agree that this step was missing from the description, and we added a paragraph for clarification.
  L184-186 It is not clear why a different correction was applied along the airborne profile corresponding to the ground data compared to the wider airborne dataset.
  We thank the reviewer for pointing this out and agree that the previous wording was unclear. We clarified that a significant local mismatch between ICECAP and both the AntGG2021 reference field, and the ground-based observation persisted after application of this regional correction. Therefore, the additional local adjustment was applied.
  L205 What observation altitude was used in the terrain model? I presume it is the same as the previously noted upward continuation altitude, but this should be stated.
  The terrain effect was calculated at the respective gravity observation elevations of each dataset. We now clarify this explicitly throughout the manuscript.
  L294 I understand you repeat the MCMC approach for each ensemble member to estimate the error, but how do you choose the ‘best’ model, or is it taken as the ensemble average of all the recovered topographies? I think this is reported on L335, but not sure.
  That is an important observation. We have clarified this and added that the representative bed model is selected as the individual realisation with the best RMS fit to the gravity data, rather than the ensemble average.
  L372-374 Suggests the gravity inversion shows topography “characterised by localised depressions and subdued ridges, which are muted or absent in current continent-scale compilations”. This is true, but are these features seen in the gravity data or artifacts of the inversion? As noted above these undulations are seen away from any data (gravity or radar). It may be that the inversion is generating realistic looking texture. This could be good as the undulations provide realistic texture, but it must be acknowledged that away from gravity and radar observations this is still a form of artificial noise.
  Thank you, again, we agree that this distinction should be made more clearly in the manuscript. We also added a paragraph to the Results section.
  L376 MAE – this is the first occurrence of this abbreviation so should be defined.
  We have added ‘mean absolute error’.
  L453-456 These lines are duplicated immediately below (L457-L460).
  Thanks for this observation. We removed the redundant lines.
  
  Citation: https://doi.org/10.5194/egusphere-2026-166-AC2

Mareen Lösing, Alan Aitken, Michael Field, Emma MacKie, and Lu Li

Viewed

Total article views: 1,978 (including HTML, PDF, and XML)

HTML	PDF	XML	Total	BibTeX	EndNote
852	1,035	91	1,978	219	268

HTML: 852
PDF: 1,035
XML: 91
Total: 1,978
BibTeX: 219
EndNote: 268

Views and downloads (calculated since 10 Feb 2026)

Month	HTML	PDF	XML	Total
Feb 2026	478	342	49	869
Mar 2026	267	534	36	837
Apr 2026	65	87	2	154
May 2026	30	65	3	98
Jun 2026	12	7	1	20

Cumulative views and downloads (calculated since 10 Feb 2026)

Month	HTML	PDF	XML	Total
Feb 2026	478	342	49	869
Mar 2026	267	534	36	837
Apr 2026	65	87	2	154
May 2026	30	65	3	98
Jun 2026	12	7	1	20

Viewed (geographical distribution)

Total article views: 1,970 (including HTML, PDF, and XML) Thereof 1,970 with geography defined and 0 with unknown origin.

Country	#	Views	%

Latest update: 17 Jun 2026

Short summary

The Denman Glacier in East Antarctica flows through one of the deepest valleys beneath the ice sheet and could contribute significantly to future sea level rise. Because this trough is extremely deep and narrow, existing radar-based maps poorly capture its shape. Using new ground-based gravity measurements with airborne data and statistical modeling, we reveal a more rugged and complex subglacial landscape than current topography products, indicating vulnerability to unstable retreat.


Total:	0
HTML:	0
PDF:	0
XML:	0