diff --git a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py index 226c2f24f..becbc8b54 100644 --- a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py +++ b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py @@ -132,7 +132,8 @@ def _remap_mali_topo(self): ice_mask, ice_density / mali_ocean_density * thickness <= sea_level - bed) grounded_mask = np.logical_and(ice_mask, np.logical_not(floating_mask)) - ocean_mask = np.logical_and(np.logical_not(ice_mask), bed < sea_level) + ocean_mask = np.logical_and(np.logical_not(grounded_mask), + bed < sea_level) lithop = ice_density * g * thickness ice_frac = xr.where(ice_mask, 1., 0.) @@ -236,17 +237,29 @@ def _combine_topo(self): 'landIceThkObserved']: ds_out[var] = ds_mali[var] - # for now, blend topography at calving fronts, but we may want a - # smoother blend in the future - alpha = ds_mali.landIceFracBilinear - ds_out['bed_elevation'] = (alpha * ds_mali.bed_elevation + - (1.0 - alpha) * ds_bedmachine.bed_elevation) + ds_out['landIceFracObserved'] = ds_mali['landIceFracConserve'] + + ds_out['landIceFloatingFracObserved'] = ( + ds_mali['landIceFracConserve'] - + ds_mali['landIceGroundedFracConserve']) for var in ['maliFracBilinear', 'maliFracConserve']: mali_field = ds_mali[var] mali_field = xr.where(mali_field.notnull(), mali_field, 0.) ds_out[var] = mali_field + # for now, blend topography at calving fronts, but we may want a + # smoother blend in the future + alpha = ds_mali.landIceFracBilinear + ds_out['bed_elevation'] = ( + alpha * ds_mali.bed_elevation + + (1.0 - alpha) * ds_bedmachine.bed_elevation) + + alpha = ds_out.maliFracConserve + ds_out['oceanFracObserved'] = ( + alpha * ds_mali.oceanFracConserve + + (1.0 - alpha) * ds_bedmachine.oceanFracObserved) + ds_out['ssh'] = ds_out.landIceDraftObserved write_netcdf(ds_out, 'topography_remapped.nc')