5

Aloha! I'm trying to plot the paths of multiple ships on a world map using R and a series of Lat/Lon points.

All is well until the ship crosses the international dateline (-180/180) when the path jumps across the map.

I tried applying the st_wrap_dateline() function in R but it looks like it only works on datasets with two points, a start & end coordinate.

My R code is below as well as a screenshot of the plot produced - any help is much appreciated!

World Map showing the errant path

# Download background Blue Marble Globe image
download.file("https://www.researchvessels.org/images/nasa_base_v2.png", "nasa_base_v2.png", mode="wb")
# Download ship lat/lon data
download.file("https://www.researchvessels.org/images/shipdata.RDATA", "shipdata.RDATA", mode="wb")
load("shipdata.RDATA", envir = parent.frame(), verbose=TRUE)

# Ensure that the required packages are installed
list.of.packages <- c("ggplot2", "sf", "dplyr", "png", "grid")
new.packages <- list.of.packages[!(list.of.packages %in% installed.packages()[,"Package"])]
if(length(new.packages)>0) {install.packages(new.packages)}


library(ggplot2)
library(sf)
library(dplyr)
library(png)
library(grid)

bluemarble_bg <- png::readPNG("nasa_base_v2.png")

xlim = c(-180.0,180.0)
ylim = c(-90.0,90.0)

shiptracks$Date.Time <- as.POSIXct(shiptracks$Date.Time, tz="GMT", origin="1970-01-01")

pal <- c("ShipA" = "#488f31", 
         "ShipB" = "#FF00FF", 
         "ShipC" = "#fff1a9",
         "ShipD" = "#f19d61", 
         "ShipE" = "#de425b")

shiptracks %>% group_by(Vessels.Name)

shiptracks <- st_as_sf(shiptracks, coords=c("Lon", "Lat")) %>% st_set_crs(4326)

shiptracks <- cbind(shiptracks, st_coordinates(shiptracks)) 

shiptracks <- shiptracks %>%
  st_sf() %>%
  st_wrap_dateline(options = c("WRAPDATELINE=YES", "DATELINEOFFSET=180"),
                   quiet = FALSE) %>%
  sf::st_sf(crs = 4326)

shipmap <- ggplot(shiptracks, aes(x=X, y=Y, group = Vessels.Name)) +
  coord_sf(xlim = xlim, ylim = ylim, expand = FALSE, 
           crs = 4326,
           datum = sf::st_crs(4326), label_graticule = waiver(),
           label_axes = waiver(), ndiscr = 100, default = FALSE,
           clip = "on") +
  annotation_custom(rasterGrob(bluemarble_bg, 
                               width = unit(1,"npc"), 
                               height = unit(1,"npc")), 
                    -Inf, Inf, -Inf, Inf) +
  geom_path(data = shiptracks, group = shiptracks$Vessels.Name, color=pal[shiptracks$Vessels.Name],
            aes(x=X, y=Y), cex=1, show.legend = TRUE) +
  geom_point(data = shiptracks, group = shiptracks$Vessels.Name, color=pal[shiptracks$Vessels.Name], aes(x=X, y=Y), cex=1, show.legend = TRUE) +
    theme(axis.line = element_blank(),
        axis.text.x = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks = element_blank(),
        axis.title.x = element_blank(),
        axis.title.y = element_blank(),
        panel.background = element_blank(),
        panel.border = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        plot.background = element_blank(),
        plot.margin = unit(c(0, 0, 0, 0), "mm"),
        panel.spacing = unit(0, "mm"),
        legend.position = c(-120, -20),
        legend.title = element_text(size = 7), 
        legend.text = element_text(size = 5), 
        legend.key.size = unit(0.01, "npc"))

shipmap
cpuguru
  • 3,763
  • 5
  • 33
  • 31
  • As more crossings happen, the visualization gets quite messy: https://www.researchvessels.org/images/arf_vis.gif – cpuguru Nov 22 '19 at 19:11

1 Answers1

5

At the point at which you run st_wrap_dateline you have a data frame of POINT geometries. Points can't cross the dateline, and so nothing changes.

You need to create a data frame of LINESTRING or MULTILINESTRING geometries by grouping and construction, so you end up with one row per ship track, and not one row per ship location with a grouping variable.

Then st_wrap_dateline will split those at the date line into MULTILINESTRING objects.

Then you can use geom_sf in your ggplot call to draw the lines, which will have been split into parts at the dateline.

Step 1, create lines from tracks:

tracks = shiptracks %>% group_by(Vessels.Name) %>% 
    summarise(do_union=FALSE) %>% st_cast("LINESTRING")

ggplot() + geom_sf(data=tracks,aes(col=Vessels.Name))

enter image description here

Now wrap it and plot:

tracks = st_wrap_dateline(tracks)
ggplot() + geom_sf(data=tracks,aes(col=Vessels.Name)) 

enter image description here

(Tip: please try and keep your code to a minimum - my first efforts at this failed because you must be using an updated version of something over what I've got, and one of those numerous options in your plots (waiver() I think) wasn't recognised. Simplifying everything sometimes also reveals where the problem is, or at least makes it quicker to find it. It also encourages people to answer. When I first saw this and its 60 lines of code I was disinclined to look further.)

Spacedman
  • 92,590
  • 12
  • 140
  • 224
  • Thanks! Any thoughts on coercing this back into the original data as additional points? The next step is for me to use the data in the Move package to animate the movements. – cpuguru Nov 25 '19 at 19:35