Skip to content

Error in st_cast.POINT(X[[i]], ...) : cannot create MULTILINESTRING from POINT when supplying named numeric vector to coord_sf() #6818

@ltyrone

Description

@ltyrone

In the following, why does hard coding ylim = work, but supplying a vector without defining the datum in coord_sf() throw an error?

I understand that the water object extent causes the issue, but it is that size to stop other map features (not shown here) from extending beyond the globe extent.

I tried y_lim <- as.numeric(...) but it did not help.

This ticket is the closest I found, but my data are projected.

Reprex:

library(ggplot2) # ggplot2_4.0.2
library(sf) # sf_1.1-0
library(rnaturalearth) # rnaturalearth_1.2.0
library(dplyr) # dplyr_1.1.4

# Create bbox
bbox <- st_as_sfc(st_bbox(c(xmin = -50,
                            ymin = 0,
                            xmax = 120,
                            ymax = 90))) |>
  st_as_sf(crs = 4326) |>
  st_segmentize(1e4)

# Centroid coords
centroid <- st_coordinates(st_centroid(bbox))

# Define PROJ4 string
prj <- paste0("+proj=ortho +lat_0=", centroid[2], " +lon_0=", centroid[1])

# Load world sf
world <- ne_countries() |>
  st_transform(prj) |>
  st_make_valid() |>
  filter(!st_is_empty(geometry))

# Create 'water' polygon
water <- st_sfc(st_point(c(35, 43.47248)), crs = 4326) |>
  st_transform(prj) |>
  st_buffer(6377341)

# Get projected bbox lat coordinates for ylim extent
y_lim <- st_bbox(st_transform(bbox, ortho_prj_35_43))[c(2,4)]

p <- ggplot() +
  geom_sf(data = water,
          fill = "#C8DAF7",
          colour ="#C8DAF7",
          linewidth = 0.25) +
  geom_sf(data = world,
          colour = "grey35",
          fill = "grey55",
          linewidth = 0.1)
  
# Both work
p + coord_sf(expand = FALSE,
           ylim = c(-4366845, 4634479)) +
  theme_void()

p + coord_sf(expand = FALSE,
             ylim = y_lim,
             datum = ortho_prj_35_43) +
  theme_void()

# Error
p + coord_sf(expand = FALSE,
             ylim = y_lim) +
  theme_void()

Error in st_cast.POINT(X[[i]], ...) : 
  cannot create MULTILINESTRING from POINT

traceback()
21: stop("cannot create MULTILINESTRING from POINT")
20: st_cast.POINT(X[[i]], ...)
19: FUN(X[[i]], ...)
18: lapply(x, st_cast, to = to)
17: st_sfc(lapply(x, st_cast, to = to), crs = st_crs(x), precision = st_precision(x))
16: st_cast.sfc(st_geometry(x), to, group_or_split = do_split)
15: st_cast(st_geometry(x), to, group_or_split = do_split)
14: st_cast.sf(df, "MULTILINESTRING")
13: st_cast(df, "MULTILINESTRING")
12: st_cast(st_cast(df, "MULTILINESTRING"), "LINESTRING", warn = FALSE)
11: sf::st_graticule(bbox, crs = params$crs, lat = breaks$y %|W|% 
        NULL, lon = breaks$x %|W|% NULL, datum = self$datum, ndiscr = self$ndiscr)
10: setup_panel_params(..., self = self)
9: (structure(function (...) 
   setup_panel_params(..., self = self), class = "ggproto_method"))(dots[[1L]][[1L]], 
       dots[[2L]][[1L]], params = list(expand = c(FALSE, FALSE, 
       FALSE, FALSE), crs = structure(list(input = "+proj=ortho +lat_0=43.4724782010486 +lon_0=35", 
           wkt = "PROJCRS[\"unknown\",\n    BASEGEOGCRS[\"unknown\",\n        DATUM[\"World Geodetic System 1984\",\n            ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n                LENGTHUNIT[\"metre\",1]],\n            ID[\"EPSG\",6326]],\n        PRIMEM[\"Greenwich\",0,\n            ANGLEUNIT[\"degree\",0.0174532925199433],\n            ID[\"EPSG\",8901]]],\n    CONVERSION[\"unknown\",\n        METHOD[\"Orthographic\",\n            ID[\"EPSG\",9840]],\n        PARAMETER[\"Latitude of natural origin\",43.4724782010486,\n            ANGLEUNIT[\"degree\",0.0174532925199433],\n            ID[\"EPSG\",8801]],\n        PARAMETER[\"Longitude of natural origin\",35,\n            ANGLEUNIT[\"degree\",0.0174532925199433],\n            ID[\"EPSG\",8802]],\n        PARAMETER[\"False easting\",0,\n            LENGTHUNIT[\"metre\",1],\n            ID[\"EPSG\",8806]],\n        PARAMETER[\"False northing\",0,\n            LENGTHUNIT[\"metre\",1],\n            ID[\"EPSG\",8807]]],\n    CS[Cartesian,2],\n        AXIS[\"(E)\",east,\n            ORDER[1],\n            LENGTHUNIT[\"metre\",1,\n                ID[\"EPSG\",9001]]],\n        AXIS[\"(N)\",north,\n            ORDER[2],\n            LENGTHUNIT[\"metre\",1,\n                ID[\"EPSG\",9001]]]]"), class = "crs")))
8: mapply(FUN = f, ..., SIMPLIFY = FALSE)
7: Map(self$coord$setup_panel_params, scales_x, scales_y, MoreArgs = list(params = self$coord_params))
6: setup_panel_params(..., self = self)
5: layout$setup_panel_params()
4: `ggplot_build.ggplot2::ggplot`(x)
3: ggplot_build(x)
2: `print.ggplot2::ggplot`(x)
1: (function (x, ...) 
   UseMethod("print"))(x)

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions