PBSmapping: importGSHHS – R documentation

Become an expert in R — Interactive courses, Cheat Sheets, certificates and more!

importGSHHS

Import Data from a GSHHS Database

Description

Import data from a GSHHS database and convert data into a PolySet with a PolyData attribute. The database was originally called ‘Global Self-consistent, Hierarchical, High-resolution Shoreline’ (GSHHS, Wessel and Smith 1996), but ‘Shoreline’ was subsequently expanded to include more ‘Geography’ (GSHHG).

Usage

importGSHHS(gshhsDB, xlim, ylim, maxLevel=4, n=0, useWest=FALSE)

Arguments

`gshhsDB`	path name to binary GSHHS database. If unspecified, looks for `gshhs_f.b` in the root of the PBSmapping library directory.
`xlim`	range of X-coordinates (for clipping). The range should be between 0 and 360, starting at the Greenwich meridian and wrapping eastward around the globe.
`ylim`	range of Y-coordinates (for clipping).
`maxLevel`	maximum level of polygons to import: 1 (land), 2 (lakes on land), 3 (islands in lakes), or 4 (ponds on islands); ignored when importing lines.
`n`	minimum number of vertices that must exist in a line/polygon in order for it to be imported.
`useWest`	logical: if `TRUE`, convert the X-coordinates (longitude) to degW (western hemisphere -180 to 0, i.e., west of the Greenwich meridian).

Details

This routine requires a binary GSHHG (Global Self-consistent, Hierarchical, High-resolution Geography) database file. The GSHHG database has been released in the public domain. At the time of writing, the most recent binary database was the archive file called gshhg-bin-2.3.7.zip.

The archive contains multiple binary files that contain geographical coordinates for shorelines (gshhs), rivers (wdb_rivers), and borders (wdb_borders). The latter two come from World DataBank II (WDBII). The five resolutions available are: full (f), high (h), intermediate (i), low (l), and coarse (c).

This routine returns a PolySet object with an associated PolyData attribute. The attribute contains four fields: (a) PID, (b) SID, (c) Level, and (d) Source. Each record corresponds to a line/polygon in the PolySet. The Level indicates the line's/polygon's level (1=land, 2=lake, 3=island, 4=pond). The Source identifies the data source (1=WVS, 0=CIA (WDBII)).

Value

A PolySet with a PolyData attribute.

Note

The function calls a C routine, also called importGSHHS, which returns a set of map coordinates that is not always predictably laid out. This issue stems from how the world is divided at the Greenwich meridian and at the International Date Line. The unpredictability occurs when user-specified X-limits span either of the longitudinal meridians – (0deg, 360deg) or (-180deg, 180deg).

This version of the R function attempts to stitch together the overlapping edges of gshhs that run from -20deg to 360deg (see example map 5 below). At present, no attempt has been made to deal with the overlap at the International Date Line where Russia overlaps the Aleutian Islands of Alaska. To some extent, the C-code can deal with this, but not in all cases.

Therefore, the user will likely experience some limitations when using importGSHHS. The solution is to import the whole dataset with this function using xlim=c(0,360), and then apply the function refocusWorld with user-desired X-limits. The Y-limits are generally not problematic unless the user wants to focus on either pole.

Author(s)

Nicholas M. Boers, Associate Professor – Computer Science
MacEwan University, Edmonton AB

Maintainer: Rowan Haigh, Program Head – Offshore Rockfish
Pacific Biological Station (PBS), Fisheries & Oceans Canada (DFO), Nanaimo BC
locus opus: Institute of Ocean Sciences (IOS), Sidney BC
Last modified Rd: 2021-01-11

References

Wessel, P., and Smith, W.H.F. (1996) A Global Self-consistent, Hierarchical, High-resolution Shoreline Database. J. Geophys. Res. 101 8741-8743.

Examples

## Not run: 
useWest=FALSE
useVers=c("2.2.0","2.2.3","2.3.0","2.3.4") # GSHHG versions
mapswitch = 5
for (i in c("land","rivers","borders"))
  if (exists(i)) eval(parse(text=paste0("rm(",i,")")))
switch( mapswitch,
 # 1. Canada------------------------------------------------
      {vN=4; useWest=T; xlim=c(-150,-50)+360;ylim=c(40,75)},
 # 2. NW Canada & America-----------------------------------
      {vN=4; useWest=T;xlim=c(-136,-100)+360;ylim=c(40,75)},
 # 3. Black Sea (user Ivailo)-------------------------------
      {vN=4; xlim=c(27.5, 34.3); ylim=c(40.9, 46.7)},
 # 4. W Europe, NW Africa (user Uli)------------------------
      {vN=4; xlim=c(-20,10); ylim=c(20,50)},
 # 5. W Europe + Iceland------------------------------------
      {vN=4; xlim=c(-25, 20); ylim=c(40, 68)},
 # 6. New Zealand-------------------------------------------
      {vN=4; xlim=c(163, 182); ylim=c(-48,-34)},
 # 7. Australia---------------------------------------------
      {vN=4; xlim=c(112,155); ylim=c(-44,-10)},
 # 8. Japan-------------------------------------------------
      {vN=4; xlim=c(127,148); ylim=c(30,47)},
 # 9. Central America---------------------------------------
      {vN=4; useWest=T; xlim=c(-95,-60)+360;ylim=c(-10,25)},
 #10. North Pacific-----------------------------------------
      {vN=4; useWest=T; xlim=c(150,220); ylim=c(45,80)},
 #11. Pacific Ocean-----------------------------------------
      {vN=4; xlim=c(112,240); ylim=c(-48,80)},
 #12. North Atlantic (world coordinates)--------------------
      {vN=4; xlim=c(285,360); ylim=c(40,68)},
 #13. North Atlantic (western hemisphere coordinates)-------
      {vN=4; xlim=c(-75,0); ylim=c(40,68)},
 #14. Atlantic Ocean----------------------------------------
      {vN=4; xlim=c(285,380); ylim=c(-50,68)},
 #15. Northern hemisphere-----------------------------------
      {vN=4; xlim=c(-180,180); ylim=c(0,85)},
 #16. Asia--------------------------------------------------
      {vN=4; xlim=c(0,180); ylim=c(0,80)},
 #17. North America-----------------------------------------
      {vN=4; xlim=c(-180,0); ylim=c(0,80)},
 #18. International date line-------------------------------
      {vN=4; xlim=c(45,315); ylim=c(0,80)},
 #19. Indian Ocean------------------------------------------
      {vN=4; xlim=c(20,130); ylim=c(-40,40)},
 #20. Moose County ("400 miles north of everywhere")--------
      {vN=4; xlim=c(272.5,280.5); ylim=c(43,47.5)}
)
db=paste0("gshhg-bin-",useVers[vN])        # database version folder
gshhg   = paste0("C:/Ruser/GSHHG/",db,"/") # directory with binary files
land    = importGSHHS(paste0(gshhg,"gshhs_i.b"),
          xlim=xlim,ylim=ylim,maxLevel=4,useWest=useWest)
rivers  = importGSHHS(paste0(gshhg,"wdb_rivers_i.b"),
          xlim=xlim,ylim=ylim,useWest=useWest)
borders = importGSHHS(paste0(gshhg,"wdb_borders_i.b"),
          xlim=xlim,ylim=ylim,useWest=useWest,maxLevel=1)
if(exists("land")){
  plotMap(land,xlim=xlim-ifelse(useWest,360,0),ylim=ylim,
    col="lemonchiffon",bg="aliceblue")
  if(!is.null(rivers)) addLines(rivers,col="blue")
  if(!is.null(borders)) addLines(borders,col="red",lwd=2)
}

## End(Not run)

PBSmapping

Mapping Fisheries Data and Spatial Analysis Tools

v2.73.0

GPL (>= 2)

Authors

Jon T. Schnute [aut], Nicholas Boers [aut], Rowan Haigh [aut, cre], Alex Couture-Beil [ctb], Denis Chabot [ctb], Chris Grandin [ctb], Angus Johnson [ctb], Paul Wessel [ctb], Franklin Antonio [ctb], Nicholas J. Lewin-Koh [ctb], Roger Bivand [ctb]

Initial release

2021-01-12