MyNameisTK
MyNameisTK

Reputation: 229

How to reduce the space between the plot and the border for geographic maps?

I am trying to plot a bathymetry map of the the northeast US using the marmap library. The following code loads the correct extent but when I plot the map I have blank space between the border and the map either at the top/bottom or left/right of the map. This also occurs when exporting the plots. If I drag the plot viewer screen size the plot adjusts and I can remove almost all of the empty space but I will be running this script in a loop so its not practical to solve this problem this way. Because of the loop I also can't hard code any dimensions into the plot because it will change for each new extent. How can I set the border of the plot to match the extent of the bathymetry?

library(marmap)
library(maps)

atl<- getNOAA.bathy(-80.93645,-41.61417,30.2 ,60.905 ,resolution=4)
blues <- colorRampPalette(c("darkblue", "cyan"))
greys <- colorRampPalette(c(grey(0.4),grey(0.99)))

plot(atl, image = TRUE, land = TRUE, n=0,
     bpal = list(c(0, max(atl), greys(100)),
                 c(min(atl), 0, blues(100))))
map(database= "state", col="black", fill=FALSE, add=TRUE)
text(x=state.center$x, y=state.center$y, state.abb, cex=0.5)

enter image description here

Upvotes: 2

Views: 1228

Answers (2)

DJack
DJack

Reputation: 4940

Here is a proposal using a workaround. The idea is to convert the bathy object into raster object and then make the plot using levelplot from rasterVisthat correctly fits the plotting area to the raster extent. Note that using raster allows having a defined pixel size and, therefore, a correct width/height ratio that you don't seem to have with marmap::plot method.

library(raster)
library(rasterVis)

r <- marmap::as.raster(atl)

state <- map('state', plot = FALSE) 
state <- data.frame(lon = state$x, lat = state$y) 

state.lab <- data.frame(lon = state.center$x, lat = state.center$y, 
                        label = state.abb)

# you can remove the color legend by adding colorkey = FALSE in levelplot() 
levelplot(r, 
          at = c(seq(min(atl), 0, length.out = 100), 
                 seq(0, max(atl), length.out = 100)[-1]),
          col.regions = c(blues(100), greys(100)), 
          margin = FALSE) +
  xyplot(lat ~ lon, state, type = 'l', 
         col = 'black') +
  xyplot(lat ~ lon, data = state.lab,
         panel = function(y, x, ...) {
         ltext(x = x, y = y, labels = state.lab$label, cex = 0.75)
         })

enter image description here

Upvotes: 2

Benoit
Benoit

Reputation: 1234

This behavior is caused by the asp argument of plot.bathy(). By default, it is fixed as asp = 1 to ensure that the scales on both axes are the same (one degree of longitude equals one degree of latitude). An unwelcome consequence of this default, is the white bands appearing either on the left/right sides of the graph, or on the top/bottom sides depending on the dimensions of your bathymetric map and the plotting device.

So I suppose you have 2 options:

  1. If you don't mind having a slightly distorted perspective, you can set asp = NA in your call to plot.bathy()
  2. If you want to have the correct aspect ratio but need to use the default size for your plotting region, then you have to download a bathymetric region that covers the whole plotting region of your active device. For instance, you could call plot.bathy() once to create a "default" plot, then, use par("usr") to determine the limits of the bathymetry needed to fill the entire plotting area. You would then download a second bathymetry with the appropriate ranges in longitude and latitude. Which is maybe not desirable.

Here is what the code would look like for the second option:

atl <- getNOAA.bathy(-80.93645, -41.61417, 30.2, 60.905, resolution = 4)
blues <- colorRampPalette(c("darkblue", "cyan"))
greys <- colorRampPalette(c(grey(0.4), grey(0.99)))

plot(atl, image = TRUE, land = TRUE, n = 0,
     bpal = list(c(0, max(atl), greys(100)),
                 c(min(atl), 0, blues(100))))

coord <- par("usr")
atl2 <- getNOAA.bathy(coord[1], coord[2], coord[3], coord[4], res = 4)
plot(atl2, image = TRUE, land = TRUE, lwd = 0.2,
     bpal = list(c(0, max(atl2), greys(100)),
                 c(min(atl2), 0, blues(100))))
map(database = "state", col = "black", fill = FALSE, add = TRUE)
text(x = state.center$x, y = state.center$y, state.abb, cex = 0.5)

enter image description here

I suppose the solution proposed by Roman Luštrik works too, but it has the inconvenience of leaving the white bands visible on both sides of the plot.

As an aside, if you have a lot of bathymetric regions to plot, you should maybe consider using the keep = TRUE argument of getNOAA.bathy() to avoid querying the NOAA servers each time you need to re-execute your code (and it is much faster to load local data than remote ones). And you could also download once and for all the global 4Go ETOPO1 and use subset.bathy() to, well, subset the bathymetry you need for each plot.

Upvotes: 2

Related Questions