Transport chapter failing on actions

[Robinlovelace]

With the following message:

Quitting from lines 376-414 (13-transport.Rmd) 
Error: arguments have different crs
In addition: Warning message:
In st_centroid.sf(zones_od) :
  st_centroid assumes attributes are constant over geometries of x

Execution halted
Error in Rscript_render(f, render_args, render_meta, add1, add2) : 
  Failed to compile 13-transport.Rmd
Calls: <Anonymous> -> render_new_session -> Rscript_render
Execution halted
Error: Process completed with exit code 1.

Source: https://github.com/Robinlovelace/geocompr/runs/8072585072?check_suite_focus=true#step:4:8607

[Robinlovelace]

Not sure if this is related but the latest data from spDataLarge is way different.

Before:

After:

And interactive map:

[Robinlovelace]

I think I've found the cause, from spDataLarge tests:

waldo::compare(bristol_stations, bristol_stations_new)
lines(attr(old$geometry, 'crs')$wkt)[6:11] vs lines(attr(new$geometry, 'crs')$wkt)[6:12]
  "        MEMBER[\"World Geodetic System 1984 (G873)\"],"
  "        MEMBER[\"World Geodetic System 1984 (G1150)\"],"
  "        MEMBER[\"World Geodetic System 1984 (G1674)\"],"
  "        MEMBER[\"World Geodetic System 1984 (G1762)\"],"
+ "        MEMBER[\"World Geodetic System 1984 (G2139)\"],"
  "        ELLIPSOID[\"WGS 84\",6378137,298.257223563,"
  "            LENGTHUNIT[\"metre\",1]],"
  "        ENSEMBLEACCURACY[2.0]],"

[Robinlovelace]

Should now be fixed upstream in spDataLarge

[Robinlovelace]

Being installed on CI, fingers crossed:

Run Rscript -e 'remotes::install_github("geocompr/geocompkg", dependencies = TRUE, force = TRUE)'

Using github PAT from envvar GITHUB_PAT
Downloading GitHub repo geocompr/geocompkg@HEAD
spDataLarge (071e59bef... -> ccc397f69...) [GitHub]
bbotk       (0.5.3        -> 0.5.4       ) [CRAN]
mlr3tuning  (0.13.1       -> 0.14.0      ) [CRAN]

Source: https://github.com/Robinlovelace/geocompr/runs/8074387299?check_suite_focus=true

[Robinlovelace]

Still failing, not sure why, and cannot reproduce locally.

[Robinlovelace]

One thing for sure: it's not due to spatial subsetting:

zone_cents_rail = zone_cents[desire_rail, ]
Error in st_geos_binop("intersects", x, y, sparse = sparse, prepared = prepared,  : 
  st_crs(x) == st_crs(y) is not TRUE

Error is different.

[Robinlovelace]

Working hypothesis: this is the culprit.

bristol_rail_points = rbind(
  st_sf(data.frame(
    Node = "Origin and destination locations",
    col = "black"
    ), geometry = zone_cents_rail$geometry),
  st_sf(data.frame(
    Node = "Public transport node",
    col = "red"
    ), geometry = bristol_stations$geometry)
)

[Robinlovelace]

This seems to be the cause:

waldo::compare(sf::st_crs(bristol_stations), sf::st_crs(bristol_zones))
lines(old$wkt) vs lines(new$wkt)
- "GEOGCRS[\"WGS 84\","
+ "GEOGCS[\"WGS 84\","
- "    ENSEMBLE[\"World Geodetic System 1984 ensemble\","
+ "    DATUM[\"WGS_1984\","
- "        MEMBER[\"World Geodetic System 1984 (Transit)\"],"
+ "        SPHEROID[\"WGS 84\",6378137,298.257223563,"
- "        MEMBER[\"World Geodetic System 1984 (G730)\"],"
+ "            AUTHORITY[\"EPSG\",\"7030\"]],"
- "        MEMBER[\"World Geodetic System 1984 (G873)\"],"
+ "        AUTHORITY[\"EPSG\",\"6326\"]],"
- "        MEMBER[\"World Geodetic System 1984 (G1150)\"],"
+ "    PRIMEM[\"Greenwich\",0,"
- "        MEMBER[\"World Geodetic System 1984 (G1674)\"],"
+ "        AUTHORITY[\"EPSG\",\"8901\"]],"
- "        MEMBER[\"World Geodetic System 1984 (G1762)\"],"
+ "    UNIT[\"degree\",0.0174532925199433,"
- "        MEMBER[\"World Geodetic System 1984 (G2139)\"],"
+ "        AUTHORITY[\"EPSG\",\"9122\"]],"
- "        ELLIPSOID[\"WGS 84\",6378137,298.257223563,"
+ "    AUTHORITY[\"EPSG\",\"4326\"]]"

[Robinlovelace]

Reprex showing it working locally:

# Aim: find out the cause of the build failing.
library(tidyverse)
remotes::install_github("nowosad/spDataLarge")
#> Skipping install of 'spDataLarge' from a github remote, the SHA1 (fe212397) has not changed since last install.
#>   Use `force = TRUE` to force installation
library(spDataLarge)
library(stplanr)
library(sf)
#> Linking to GEOS 3.10.2, GDAL 3.4.3, PROJ 8.2.0; sf_use_s2() is TRUE
library(tmap)

od_inter = filter(bristol_od, o != d)
desire_lines = od2line(od_inter, bristol_zones)
#> Creating centroids representing desire line start and end points.
desire_rail = top_n(desire_lines, n = 3, wt = train)
desire_rail = line_via(desire_rail, bristol_stations)

# Relevant lines failing:
zone_cents = st_centroid(bristol_zones)
#> Warning in st_centroid.sf(bristol_zones): st_centroid assumes attributes are
#> constant over geometries of x
zone_cents_rail = zone_cents[desire_rail, ]
bb = tmaptools::bb(desire_rail, ext = 1.1)
desire_rail_plot = rbind(
  st_sf(data.frame(Geometry = "Desire line (original)"), geometry = desire_rail$geometry),
  st_sf(data.frame(Geometry = "Leg 1 (origin to station)"), geometry = desire_rail$leg_orig),
  st_sf(data.frame(Geometry = "Leg 2 (station to station)"), geometry = desire_rail$leg_via),
  st_sf(data.frame(Geometry = "Leg 3 (station to destination)"), geometry = desire_rail$leg_dest)
) 
desire_rail_plot = desire_rail_plot |> 
  mutate(lty = case_when(Geometry == "Desire line (original)" ~ 2, TRUE ~ 1)) |> 
  mutate(size = case_when(Geometry == "Desire line (original)" ~ 1, TRUE ~ 2))
bristol_rail_points = rbind(
  st_sf(data.frame(
    Node = "Origin and destination locations",
    col = "black"
  ), geometry = zone_cents_rail$geometry),
  st_sf(data.frame(
    Node = "Public transport node",
    col = "red"
  ), geometry = bristol_stations$geometry)
)
tm_shape(desire_rail_plot, bbox = bb) +
  # legend incorrect in tmap v3.0.0: https://github.com/r-tmap/tmap/issues/672
  # tm_lines(col = "Geometry", palette = "Set2", lty = desire_rail_plot$lty) +
  tm_lines(col = "Geometry", palette = "Set2", lwd = "size", scale = 3, legend.lwd.show = FALSE) +
  tm_shape(bristol_rail_points) +
  # Try with different alpha values:
  # tm_dots(col = "col", size = 0.3, alpha = 0.2)
  tm_dots(col = "col", size = 0.3) +
  tm_scale_bar()

^{Created on 2022-08-29 with reprex v2.0.2}

[Robinlovelace]

Another reprex:

library(spDataLarge)
u = "https://github.com/Nowosad/spDataLarge/raw/master/data/bristol_stations.rda"
f = basename(u)
download.file(u, f)
load(f)
waldo::compare(sf::st_crs(bristol_zones), sf::st_crs(bristol_stations))
#> lines(old$wkt) vs lines(new$wkt)
#> - "GEOGCS[\"WGS 84\","
#> + "GEOGCRS[\"WGS 84\","
#> - "    DATUM[\"WGS_1984\","
#> + "    ENSEMBLE[\"World Geodetic System 1984 ensemble\","
#> - "        SPHEROID[\"WGS 84\",6378137,298.257223563,"
#> + "        MEMBER[\"World Geodetic System 1984 (Transit)\"],"
#> - "            AUTHORITY[\"EPSG\",\"7030\"]],"
#> + "        MEMBER[\"World Geodetic System 1984 (G730)\"],"
#> - "        AUTHORITY[\"EPSG\",\"6326\"]],"
#> + "        MEMBER[\"World Geodetic System 1984 (G873)\"],"
#> - "    PRIMEM[\"Greenwich\",0,"
#> + "        MEMBER[\"World Geodetic System 1984 (G1150)\"],"
#> - "        AUTHORITY[\"EPSG\",\"8901\"]],"
#> + "        MEMBER[\"World Geodetic System 1984 (G1674)\"],"
#> - "    UNIT[\"degree\",0.0174532925199433,"
#> + "        MEMBER[\"World Geodetic System 1984 (G1762)\"],"
#> - "        AUTHORITY[\"EPSG\",\"9122\"]],"
#> + "        MEMBER[\"World Geodetic System 1984 (G2139)\"],"
#> - "    AUTHORITY[\"EPSG\",\"4326\"]]"
#> + "        ELLIPSOID[\"WGS 84\",6378137,298.257223563,"
#> and 16 more ...
bristol_zones = bristol_zones |> sf::st_set_crs("EPSG:4326")
waldo::compare(sf::st_crs(bristol_zones), sf::st_crs(bristol_stations))
#> ✔ No differences

^{Created on 2022-08-29 with reprex v2.0.2}

[Robinlovelace]

Latest error:

Quitting from lines 526-541 (13-transport.Rmd) 
Error: arguments have different crs
In addition: Warning messages:
1: In st_centroid.sf(zones_od) :
  st_centroid assumes attributes are constant over geometries of x
2: In CPL_transform(x, crs, aoi, pipeline, reverse, desired_accuracy,  :
  GDAL Error 1: PROJ: proj_as_wkt: DatumEnsemble can only be exported to WKT2:2019
3: In CPL_transform(x, crs, aoi, pipeline, reverse, desired_accuracy,  :
  GDAL Error 1: PROJ: proj_as_wkt: DatumEnsemble can only be exported to WKT2:2019

Source: https://github.com/Robinlovelace/geocompr/runs/8081214010?check_suite_focus=true

I'm confused because cannot reproduce locally...

[Robinlovelace]

That was fixed, not its...

Quitting from lines 691-695 (13-transport.Rmd) 
Error in geos_op2_geom("difference", x, y, ...) : 
  st_crs(x) == st_crs(y) is not TRUE
Calls: local ... st_difference.sf -> geos_op2_df -> geos_op2_geom -> stopifnot
In addition: Warning messages:
1: In st_centroid.sf(zones_od) :
  st_centroid assumes attributes are constant over geometries of x
2: In CPL_transform(x, crs, aoi, pipeline, reverse, desired_accuracy,  :
  GDAL Error 1: PROJ: proj_as_wkt: DatumEnsemble can only be exported to WKT2:2019
3: In CPL_transform(x, crs, aoi, pipeline, reverse, desired_accuracy,  :
  GDAL Error 1: PROJ: proj_as_wkt: DatumEnsemble can only be exported to WKT2:2019
4: st_crs<- : replacing crs does not reproject data; use st_transform for that 
5: In CPL_transform(x, crs, aoi, pipeline, reverse, desired_accuracy,  :
  GDAL Error 1: PROJ: proj_as_wkt: DatumEnsemble can only be exported to WKT2:2019
Execution halted

[Robinlovelace]

It's strange because when I test the relevant objects:

> waldo::compare(
+   sf::st_crs(route_network_scenario),
+   sf::st_crs(existing_cycleways_buffer)
+ )
✔ No differences

There is no difference. Relevant lines that are failing, any ideas @Nowosad or @jannes-m ? I'm not sure but wondering if the CRS of some outputs depend on PROJ/GDAL versions.

https://github.com/Robinlovelace/geocompr/blob/7552bc6858529e668818ba2007547ec00d607ddf/13-transport.Rmd#L691-L694

[Robinlovelace]

Finally the build works! But with many work-arounds that I would like to debug and remove, so leaving open for now. As I say any thoughts on this welcome: I am bamboozled by what's happened here.

[Nowosad]

? I'm not sure but wondering if the CRS of some outputs depend on PROJ/GDAL versions.

Yes. It definitely can.