Preprocessing Shorelines and Waterbodies from GSHHS

See also:


 * Preprocessing Shorelines and Waterbodies from SWBD

Reading GSHHS data
This does not include the recently added river and border lines from CIA WDBII.

Create a table with fields of interest:

The last two commands add the geometry column and add an index on the geometry field.

GSHHS is shipped with a short C program gshhs suitable to convert the binary packed data into text format. We do some slight modifications in order to get a format that can be read by PostgreSQL. Run

Making tiles for the polar regions
Since SWBD data only cover the range between the 60th degrees of latitude, we need to complete it from GSHHS data.

Preparing the polygones
There are still some inconsistencies in polygones that are regarded as invalid by PostGIS. At least as much as polar regions are affected, it is sufficient to run a Simplify with very low tolerance to fix the polygones:

Polygons that are completly north or south of the 60th degree can be copied directly into a new table:

Now we have to intersect the polar region with all other polygones that have an overlap with polar regions. This seems to concern only the northern hemisphere. Polygons on the southern hemisphere seem not to cross the 60th degree.

Finally, we make a cut through the antarctic continent from the coast to the south pole.

Cropping into tiles
Then we crop all GSHHS polygones into 1°x1° tiles with a 0.1° buffer. The tiles are stored in a new table defined as

A short PERL script mkGshhsTileswill help to write the SQL commands for filling this table.

For the southern polar, simply flip the sign of the values for $s in the for-loop. Then run

This will run for some houres.

Converting the tiles into SWBD format
The tiles are not yet in a SWBD compatible format. GSHHS uses a hierarchical level scheme that indicates whether a polygon describes land (level 1), a lake (level 2), an island in a lake (level 3) or a pond on an island in a lake (level 4). This avoids polygons with inner rings but requires huge polygons for continets. SWBD uses polygons with inner rings and the area of the polygon are always water bodies or ocean. Each polygon has a table field named facc_code that indicates ocean, lake or river.

Ocean tiles
Let's start with ocean tiles. This is the inverse of the union of all GSHHS land tiles with respect to the tile's total area. Again, a short PERL script helps to write the SQL commands.

Again, this will run for some houres. Finally, DELETE all empty tiles.

Lake tiles
Lake polygones must not be inverted, but all island in the lakes must be used as inner rings. For lakes with islands, we run

We are careful and invent the new facc_code "BH081" (SWBD lakes have "BH080") for being able to find our lakes easily. This may be helpful for checking the result. We also may want to crop the tile's buffer from 0.1 degrees to something less.

Lakes without islands seem to be easier, but in fact, they are even more complicated since we have to find polygons that do not contain islands. Negative conditions always tend to be a mess in SQL. Fortunately, the number of lakes and islands is not too high, so we can do it like this:

Again, we have used an own facc_code. In this nested query, we first search for all lake IDs that contain islands. This is just the same as above. Then we exclude this IDs and insert the remaining lake polygons into the swbd table.

Ponds on islands can simply be copied.

Polygon orientation
The last step is to make sure all polygons are oriented left handed, i.e. the area of the polygon must be on the left hand side when walking along its rings.

When everthing looks all right, we should update the facc_code