Interesting timing of this post! I was thinking today what would be the simplest System Dynamics model that would simulate melting the ice caps of Mars using nuclear bombs. (How many, for how long, to what effect, etc.)

@Vitaliy, this is actually nice to quickly draw different depths of the North Sea. In Aberdeen we are proceeding with what they call "energy transition". For example when they are looking for $CO_2$ storage sites there are many requirements, e.g. regarding the depth of the sea floor etc. With your code you can mark different depth ranges, e.g. in the North Sea from 80-100 m like so

(* =====North Sea+NE Atlantic:relief land,blue sea,80\[Dash]\
100m in red=====*)(*Region as {{latMin,latMax},{lonMin,lonMax}}.\
Extended N/NW so a genuine 3\[Dash]4 km depth band exists.*)
region = {{50, 70}, {-14, 16}};

(*Depth band to flag,in metres (negative=below sea surface).*)
{dDeep, dShallow} = {-100, -80};   (*sea floor 80-100m down*)

(*Land ramp:green lowland->tan->brown->white peaks (metres).*)
land[e_] := 
  Blend[{RGBColor[0.24, 0.45, 0.22], RGBColor[0.45, 0.55, 0.27], 
    RGBColor[0.75, 0.68, 0.42], RGBColor[0.52, 0.39, 0.26], White}, 
   Rescale[Clip[e, {0, 2500}], {0, 2500}]];

(*Sea ramp:pale shelf->deep blue (metres).*)
sea[e_] := 
  Blend[{RGBColor[0.78, 0.88, 0.96], RGBColor[0.30, 0.55, 0.80], 
    RGBColor[0.08, 0.25, 0.55], RGBColor[0.02, 0.08, 0.25]}, 
   Rescale[Clip[e, {-5000, 0}], {0, -5000}]];

(*Master colour function.\
Argument is elevation in METRES because ColorFunctionScaling->\
False below.\
QuantityMagnitude guards against the value arriving as a Quantity in s\
ome versions.*)
elevColor = 
  Function[e0, 
   With[{e = QuantityMagnitude[e0]}, 
    Which[e >= 0, land[e],(*land*)dDeep <= e <= dShallow, 
     Red,(*3\[Dash]4 km deep sea*)True, sea[e]     (*all other sea*)]]];

GeoGraphics[{}, GeoRange -> region, GeoProjection -> "Mercator", 
 GeoBackground -> 
  GeoStyling["ReliefMap", ColorFunction -> elevColor, 
   ColorFunctionScaling -> False], GeoGridLines -> Automatic, 
 ImageSize -> 850]

If you want to shade the different depths, that looks like this:

    (* =====North Sea+NE Atlantic:relief land,blue sea,80\[Dash]\
    100m in red=====*)(*Region as {{latMin,latMax},{lonMin,lonMax}}.\
    Extended N/NW so a genuine 3\[Dash]4 km depth band exists.*)
    region = {{50, 70}, {-14, 16}};

    (*Depth band to flag,in metres (negative=below sea surface).*)
    {dDeep, dShallow} = {-100, -80};   (*sea floor 80-100m down*)

    (*Land ramp:green lowland->tan->brown->white peaks (metres).*)
    land[e_] := 
      Blend[{RGBColor[0.24, 0.45, 0.22], RGBColor[0.45, 0.55, 0.27], 
        RGBColor[0.75, 0.68, 0.42], RGBColor[0.52, 0.39, 0.26], White}, 
       Rescale[Clip[e, {0, 2500}], {0, 2500}]];

    (*Sea ramp:pale shelf->deep blue (metres).*)
    sea[e_] := 
      Blend[{RGBColor[0.78, 0.88, 0.96], RGBColor[0.30, 0.55, 0.80], 
        RGBColor[0.08, 0.25, 0.55], RGBColor[0.02, 0.08, 0.25]}, 
       Rescale[Clip[e, {-5000, 0}], {0, -5000}]];

    (*Master colour function.\
    Argument is elevation in METRES because ColorFunctionScaling->\
    False below.\
    QuantityMagnitude guards against the value arriving as a Quantity in s\
    ome versions.*)
    elevColor = 
      Function[e0, 
       With[{e = QuantityMagnitude[e0]}, 
        Which[e >= 0, land[e],(*land*)dDeep <= e <= dShallow, 
         Red,(*3\[Dash]4 km deep sea*)True, sea[e]     (*all other sea*)]]];

    GeoGraphics[{}, GeoRange -> region, GeoProjection -> "Mercator", 
     GeoBackground -> 
      GeoStyling["ReliefMap", ColorFunction -> elevColor, 
       ColorFunctionScaling -> False], GeoGridLines -> Automatic, 
     ImageSize -> 850]


![enter image description here][2]

It turns out that the North Sea is really shallow....

(*North Sea shelf:sea floor 0 to-100 m,coloured by depth*)
region = {{50, 62}, {-5, 10}};         (*North Sea*)
dMax = 100;                           (*colour scale spans 0 to 100 m*)
step = 10;                            \
(*band thickness,m;set 0 for smooth*)
seaScheme = 
  ColorData[
   "TemperatureMap"];   \
(*or "Rainbow","DeepSeaColors","SouthwestColors"*)

land[e_] := 
  Blend[{RGBColor[0.24, 0.45, 0.22], RGBColor[0.45, 0.55, 0.27], 
    RGBColor[0.75, 0.68, 0.42], RGBColor[0.52, 0.39, 0.26], White}, 
   Rescale[Clip[e, {0, 2500}], {0, 2500}]];

depthColor = 
  Function[e0, 
   With[{e = QuantityMagnitude[e0]}, 
    If[e >= 0, land[e], 
     seaScheme[
      Clip[(If[step > 0, step*Floor[(-e)/step], -e])/dMax, {0, 1}]]]]];

Legended[
 GeoGraphics[{}, GeoRange -> region, GeoProjection -> "Mercator", 
  GeoBackground -> 
   GeoStyling["ReliefMap", ColorFunction -> depthColor, 
    ColorFunctionScaling -> False], GeoGridLines -> Automatic, 
  ImageSize -> 750, 
  PlotLabel -> Style["North Sea shelf \[LongDash] 0 to 100 m", 14]], 
 BarLegend[{seaScheme, {0, dMax}}, LegendLabel -> "depth (m)", 
  LegendMarkerSize -> 320]]

  [1]: https://community.wolfram.com//c/portal/getImageAttachment?filename=NorthSea80-100mdepth.jpg&userId=48754

Hi @Vitaliy Kaurov ,

this is great! I had been thinking of doing this for some other project, but never got around to do it. It is really nice that this works with so little code.

I have tried this out for our home planet for three different "filling ranges".

1. Earth from 0 to 1

Earth is filled from level 0 to 1, i.e from no water up to everything is filled. You notice that only in the very last frame(s) we reach from the current state to a completely watery planet.

1. Earth with a zoom from 0.99 to 1

This is just a look at the later filling stages.

1. Earth filled from 0.997 to 1

This looks more like a state from a bit less water than now to "all filled with water".

By the way, if you do that for the moon you see the effect of rather deep impact craters.

I guess that bodies with little to no atmosphere and smaller bodies will have a more interesting pattern.

So if you start to terra-form Mercury (and I think it might be difficult for it to retain any water) we get this here:

Thank you for this post and giving me the opportunity to play a bit with your code.

Thanks, Marco

PS: 71% surface water level should not be misunderstood as 71% of the surface. You describe it well as "fill the bottom 26% (or "71%" of the planet’s total height", or height variation.

Attachments:

earthZoom1.gif

earthfromsuperzoom.gif

earthfromzero.gif

mercury.gif

moon.gif

Hi @Vitaliy Kaurov

I used your exact code:

Monitor[frames4 = 
   Table[elSamp2 = 
     Flatten@QuantityMagnitude@
       GeoElevationData[GeoProjection -> "CylindricalEqualArea", 
        GeoZoomLevel -> 1, GeoRange -> "World", GeoModel -> "Moon"]; 
    tHeight = Rescale[Quantile[elSamp, k], MinMax[elSamp]]; 
    GeoGraphics[GeoModel -> "Moon", GeoRange -> "World", 
     GeoProjection -> "VanDerGrinten", 
     GeoBackground -> 
      GeoStyling["ReliefMap", 
       ColorFunction -> (If[# < tHeight, StandardBlue, 
           StandardOrange] &)]], {k, 0, 1, 0.01}];, k]

and then

FrameListVideo[frames4, FrameRate -> 10]

So there is no change whatsoever.

Great post!