MATLAB: Polygon areas output from the gshhs matlab function

gshhggshhspolygon areas

I'm using the gshhs function in the mapping toolbox to plot some coastlines, using the latest version of the GSHHG database maintained by Paul Wessel:
 s = gshhs('gshhs_i.b',[32 33], [-65 -64]);
This call to gshhs reads in the coastline of Bermuda. s is a 5×1 structure array that consists of fields describing five polygons: one that denotes the North American mainland and four that denote Bermuda. I assume North America is included because Bermuda is grouped with North America in the database.
I plotted the resulting polygons and verified that the last four make up Bermuda and the first makes up the North American mainland, so the lat/lons included in the structure array are what I think they are. But when I print the "Area" field included in the structure array to query the area of each polygon (in km^2, as described in the gshhs function documentation), I get the following results:
area={s.Area}
area = 
    [2.0153e+07]    [5.5061e+07]    [8.8125e+07]    [5.0053e+07]    [6.1060e+06]
The actual total area of all islands in Bermuda is about 50 km^2 (which should roughly be the sum of polygons 2-5), so I must be making a bonehead move somewhere. I would much appreciate it if someone pointed out how to get the correct polygon areas. Thank you.

Best Answer

After some further experimentation, I was able to figure out what was going on, and I thought I'd share for anyone who encounters this issue. In brief, if you want the polygon areas returned from the gshhs function to be correct, you'll need to use version 8 (aka version 2.1.1) or earlier of the GSHHG database. You can find those earlier versions here.
If you try to read in a new version of the database with gshhs, you'll receive a warning message:
 s = gshhs('gshhs_i.b',[32 33], [-65 -64]);
 Warning: The version number obtained from the file, 'gshhs_i.b', is 15. This version of the function GSHHS
 is qualified to read data up through GSHHS version 9. Attempting to continue to read the data from this version of the GSHHS dataset ...
I interpreted this warning message to mean that the function was compatible with GSHHG version 9 or earlier, and when I used version 9, the message disappeared, which is why I didn't mention it in my original post. As I reported, the polygon areas were still incorrect with v9. It turns out that version 8 is the most recent usable version if you want correct polygon areas.
Related Question