# Generate the complement of the Earth's land mass polygons

I would like to get polygon covering of all water area of earth. I though of doing by inverting

Polygon[EntityClass["Country", "LandMasses"]]


How do i invert this polygon, such that it covers all the oceans and seas?

Specific use case example:

GeoListDensityPlot[data_] := Module[
{colorRange, colorRegion, color},
colorRange = MinMax[data];
color = ColorData["TemperatureMap"];
colorRegion[pos_ -> value_] := {
color[(value - colorRange[[1]])/(colorRange[[2]] -
colorRange[[1]])], area[pos]};
GeoGraphics[{
colorRegion /@ Normal[data],
EdgeForm[Black], GeoStyling["ReliefMap"],
Polygon[EntityClass["Ocean", "Oceans"]]
} , GeoBackground -> GeoStyling["Coastlines"] ]
]


where i need to include seas as well

• How are you intending to use the resulting polygon? Do you just need to for mapping or for some sort of computation? – Edmund Aug 13 '17 at 19:14
• added a specific example, where the seas are missing in the plot – Karolis Aug 13 '17 at 19:36
• Would using GeoBackground -> GeoStyling[LightBlue] suffice? – Chip Hurst Aug 14 '17 at 14:22
• @ChipHurst data and the notebook is available on github.com/kmisiunas/data-sunshine-vs-population – Karolis Aug 14 '17 at 14:41
• [This comment is unrelated to inverting a Polygon, but my guess is you're wanting to just make a nice display. Correct me if I'm wrong...]. For the purposes of display, I don't see why you can't just set the GeoBackground to a color. Also it seems you have a regular grid of data points. If you're wanting to color the background, would making an Image be a better choice? – Chip Hurst Aug 14 '17 at 15:00

Polygon itself doesn't allow for holes and inverting a polygon will usually result in holes.

Two ways around this are with FilledCurve or BoundaryMeshRegion.

## FilledCurve

This requires some manual labor. Here's the area I'd like to invert:

polygon = CountryData["World", "SchematicPolygon"];

GeoGraphics[polygon]


world = Line[GeoPosition[{{-90, -180}, {-90, 180}, {90, 180}, {90, -180}}]];
holes = List /@ Line /@ GeoPosition /@ pp[[1, 1]];
reg = FilledCurve[Prepend[holes, {world}]];

GeoGraphics[{Red, reg}, GeoRange -> "World"]


## BoundaryMeshRegion

We can start out with a filled in earth $\left([-180, 180] \times [-90, 90]\right)$ and use RegionDifference to subtract away the polygons. I'm converting everything to mesh regions because RegionDifference usually will respect open/closed-ness of a point set, but ignores this for meshes.

For whatever reason, BoundaryDiscretizeGraphics is slow for large polygons. But Nonetheless

world = BoundaryDiscretizeGraphics[Rectangle[{-180, -90}, {180, 90}]];

holes = BoundaryDiscretizeGraphics /@ EntityValue[EntityClass["GeographicRegion", "Continents"], "Polygon"];

diff = Fold[RegionDifference, world, holes]


This is a pretty large mesh, which might be a reason why it was so slow:

MeshCellCount[diff]

{241629, 241629}


Let's drop that stray horizontal line at the bottom:

First[MaximalBy[ConnectedMeshComponents[diff], Area]]


• It's a pity that the generated mesh is huge, it looks fine otherwise. – J. M. will be back soon Aug 15 '17 at 14:53

The answer is almost there in the OP. The "Ocean" Entity is used there but only the "Oceans" EntityClass is used. There are more classes of water bodies needed other than "Oceans".

If you are not particular about overlapping polygons you may use EntityClass["Ocean", {"HasPolygon" -> True}].

GeoGraphics[{
GeoStyling[],
EdgeForm[Black],
Opacity[.5, Red],
EntityValue[EntityClass["Ocean", {"HasPolygon" -> True}], "Polygon"]
}]


As mentioned you do get some overlaps as the both the North and South polygons of the Atlantic and Pacific Oceans cover the same area as the larger (not split) Atlantic and Pacific Ocean polygons.

The can be avoided by not including the "SevenSeas" EntityClass. Unfortunately, excluding this class also excludes the oceans. Therefore, I included all classes except for the "SevenSeas" class. One quirk of this is that the "NorthwestPassage" is not in any class so it has to be included separately.

GeoGraphics[{
GeoStyling[],
EdgeForm[Black],
Opacity[.5, Red],
OceanData["NorthwestPassages", "Polygon"],
EntityValue[
EntityClass["Ocean",
{"HasPolygon" -> True,
"Memberships" ->
ContainsAny[{"Oceans", "Seas", "Bays", "Basins", "Straits", "Channels", "Gulfs"}]}],
"Polygon"]
}]


There is still some overlap in the Mediterranean Sea where it is overlapped by the polygons for is Western and Eastern basins. It should not be too difficult for you to remove these to entities before graphing if this overlap is of concern.

Hope this helps.

• This is a great answer. Thank you. I am still curious if anyone can produce an inverse from landmass polygon – Karolis Aug 13 '17 at 21:27

Just for fun. Here's a way to visualize your data and remove the small islands using image processing.

colorRange = MinMax[sunshine];
color = ColorData["Rainbow"];
cf[value_] := color[(value - colorRange[[1]])/(colorRange[[2]] - colorRange[[1]])]

im = Image[Reverse[
Partition[If[#2, cf[#1], Black] & @@@ Values[Merge[{sunshine, landQ}, Join]], 360]]];


GeoGraphics[{EdgeForm[Black], FaceForm[],