4
$\begingroup$

What is the easiest way to obtain MeshRegion objects from CountryData that I can use in geographic computations?

For example, this does not work for me:

reg = DiscretizeGraphics @ CountryData[Entity["Country", "UnitedStates"], "Polygon"];
RegionBounds[reg]

(* {{-2.75903*10^6, 1.7744*10^6}, {-5.70806*10^6, -3.48795*10^6}, {2.69158*10^6, 4.8185*10^6}} *)

I want to get mesh regions like the hardcoded USA region in the "Applications" section in RegionMember's function page. (I.e. such mesh objects for other countries.)

With the referenced hardcoded region object, reg, we get:

RegionBounds[reg]
(* {{-124.733, -66.9498}, {25.1246, 49.3845}} *)
Improve this question
$\endgroup$

2 Answers 2

Reset to default
7
$\begingroup$

Take any geo entity with polygon information and get its geo polygon:

geopol = Entity["Country", "Greenland"]["Polygon"]

The result is a Polygon[GeoPosition[...]] object that represents the actual (non-projected) shape of the region, and can be used in any of the Geo functions, in particular in GeoGraphics.

A 2D MeshRegion object with no geo information requires a choice of projection. Such object will be usable only in maps of the chosen projection. We can project the geo polygon (say to the Mercator projection) with

mercpol = GeoGridPosition[geopol, "Mercator"]

and then you can construct a 2D MeshRegion object in several ways, say

mercmesh = DiscretizeGraphics[mercpol]

Then this fits well in a Mercator map, but not for example in an equirectangular map:

{GeoGraphics[{Red, mercmesh}, GeoRange -> "World", GeoProjection -> "Mercator"],
 GeoGraphics[{Red, mercmesh}, GeoRange -> "World", GeoProjection -> "Equirectangular"]}

enter image description here

However, the original geo polygon will be correctly placed in any projection:

{GeoGraphics[{GeoStyling[Red], geopol}, GeoRange -> "World", GeoProjection -> "Mercator"],
 GeoGraphics[{GeoStyling[Red], geopol}, GeoRange -> "World", GeoProjection -> "Equirectangular"]}

enter image description here

Improve this answer
$\endgroup$
1
  • $\begingroup$ Thank you, Jose -- very insightful! $\endgroup$
    – Anton Antonov
    May 2, 2022 at 11:51
5
$\begingroup$

The key is using the cartographic projection "Equirectangular".

r4 = CountryData["UnitedStates", {"Polygon", "Equirectangular"}]
Head[r4]
r5 = MeshRegion[r4]
Region[r5, AspectRatio -> Automatic, Frame -> True]

enter image description here

RegionBounds[r5]

{{-124.733, -66.9498}, {25.1246, 49.3845}}

Improve this answer
$\endgroup$
2
  • $\begingroup$ Thank you for your prompt response! I accepted your answer and then switched the acceptance to Jose's answer -- that one explains the motivation behind the observed behavior and it is more "universal." $\endgroup$
    – Anton Antonov
    May 2, 2022 at 11:50
  • 1
    $\begingroup$ I definitely learned from it. RegionBounds[mercmesh] produces a slightly different result for "Mercator" and "Equirectangular" so I still need to learn more about how Mma handles this. $\endgroup$
    – Syed
    May 2, 2022 at 12:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.