Extraction by raster suddenly not working

[cjcarlson]

Hey Jeff!

I'm getting an error I've never gotten before using code that extracts elevation with a raster called 'rodents':

Error in (function (classes, fdef, mtable)  : 
  unable to find an inherited method for function ‘spTransform’ for signature ‘"RasterLayer", "CRS"’

I just upgraded to R 4.0 and I wonder if either the rgdal or sp versions of spTransform aren't handling things right anymore in 4? It looks like this error might be showing up in other people's auto-compiling code around the web, but I still don't entirely get what's going on:

https://rstudio-pubs-static.s3.amazonaws.com/328065_b389bbaca9ec4068b57c489d215d89e5.html

Any clues? Thanks!

[jhollist]

Sounds like that has to do with all the changes going on with the handling of spatial reference data. 'rgdal', 'sp', and 'sf' have all seen recent changes and raster may have as well. It all stems from changes in PROJ. They are good changes, but promise to wreak havoc for a while! The gist is that proj4strings are going away in favor of well-known text (WKT) representation. I think the way EPSG codes are handled has also changed. More details at https://www.r-spatial.org/r/2020/03/17/wkt.html

Can you share the code that is throwing the error? I'll see if I see anything obvious and if it is an elevatr issue, I can update on my end.

And glad to hear you are using elevatr!

[ACheysson]

Dear @jhollist,

Thanks very much for the package, it looks awesome. However I'm blocked with the same problem as @cjcarlson. Here is my code

# Load packages
library(RStoolbox)
library(raster)
library(rgdal)
library(elevatr)

# Paths to landsat-5 bands metadata 
filemeta <- 'LT05_L1TP_186031_19840910_20170220_01_T1_MTL.txt'

# Import meta-data and bands based on MTL file
metaData <- readMeta(filemeta)

# Stack the metadata bands
lsMeta <- stackMeta(metaData)

# Correct DN to at-surface-reflecatance with DOS (Chavez decay model)
l8_boa_ref <- radCor(lsMeta, metaData, method = "dos")

#Save the raster
writeRaster(l8_boa_ref,"myStack.tif", format="GTiff") 

prj_dd <- "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"

elevation_map<- get_elev_raster(l8_boa_ref@layers[[1]],prj_dd, z = 10)

I'm getting then the same bug as he does. I was just wondering, could you indicate the previous versions of sp/rgdal/sf that had elevatr working properly?

Best

[jsta]

This very simple example works for me with the absolute latest package versions:

# library(stars)
library(raster)
library(elevatr)

r <- raster(system.file("tif/L7_ETMs.tif", package="stars"))
elevation_map <- get_elev_raster(r, proj4string(r), z = 10)

session info

``` dplyr::select( dplyr::filter(as.data.frame(session_info()$packages), package %in% c("sf", "sp", "rgdal", "elevatr", "raster", "stars")), package, loadedversion, date, source) package loadedversion date source elevatr elevatr 0.3.3.9999 2021-01-07 Github (jhollist/elevatr@c772fd6) raster raster 3.4-5 2020-11-14 CRAN (R 4.0.3) rgdal rgdal 1.5-19 2021-01-05 CRAN (R 4.0.3) sf sf 0.9-7 2021-01-06 CRAN (R 4.0.3) sp sp 1.4-5 2020-12-21 Github (edzer/sp@94a5eb5) sf::sf_extSoftVersion() GEOS GDAL proj.4 GDAL_with_GEOS USE_PROJ_H "3.8.0" "3.0.4" "6.3.1" "true" "true" ```

[ACheysson]

That's great thank you very much ! just had to update the packages correctly.

[jhollist]

Thanks for testing this, @jsta ! I can cofirm that latest versions work. I am working on a solution to push to CRAN. Not much different than what is currently on GitHub. So remotes::install_github("jhollist/elevatr") should fix for now. Will update this issue once I have the new version on CRAN. Thank you all!!

[jhollist]

And it is now on CRAN!

https://cran.r-project.org/web/packages/elevatr/index.html

I'll leave this open for a week or so. Please test and confirm that the new version works for your particular uses. Thanks!

[cjcarlson]

It works! You're a hero. Thanks man!

[cjcarlson]

Ah I spoke too soon I'm sorry to say. Here's a reproducible example:

library(elevatr)
library(raster)

sessionInfo()

rodents <- raster('RodentRichness.tif')
rodents <- trim(rodents)
elev <- elevatr::get_elev_raster(rodents,z=6)

This works on my laptop but not the RStudio server I'm working on, which has switched to a new error:

> elev <- elevatr::get_elev_raster(rodents,z=6)
Error in loc_check(locations, prj) : Please supply a valid WKT string.
In addition: Warning messages:
1: In sp::wkt(locations) : no wkt comment
2: In sp::wkt(locations) : no wkt comment
3: In sp::wkt(locations) : no wkt comment

sessionInfo() on server:

> sessionInfo()
R version 4.0.2 (2020-06-22)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 16.04.7 LTS

Matrix products: default
BLAS:   /usr/lib/libblas/libblas.so.3.6.0
LAPACK: /usr/lib/lapack/liblapack.so.3.6.0

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C               LC_TIME=en_US.UTF-8       
 [4] LC_COLLATE=en_US.UTF-8     LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                  LC_ADDRESS=C              
[10] LC_TELEPHONE=C             LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] raster_3.4-5  sp_1.4-4      elevatr_0.3.3

loaded via a namespace (and not attached):
[1] compiler_4.0.2   rgdal_1.5-19     tools_4.0.2      Rcpp_1.0.5       codetools_0.2-18 grid_4.0.2      
[7] lattice_0.20-41 

sessionInfo() on laptop:

> sessionInfo()
R version 4.0.3 (2020-10-10)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 18363)

Matrix products: default

locale:
[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United States.1252   
[3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C                          
[5] LC_TIME=English_United States.1252    

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] raster_3.4-5  sp_1.4-4      elevatr_0.3.3

loaded via a namespace (and not attached):
[1] compiler_4.0.3   tools_4.0.3      Rcpp_1.0.5       codetools_0.2-16
[5] grid_4.0.3       lattice_0.20-41 

rodent file to play with:

Layers.zip

[jhollist]

Thanks for the example. I'll take a look, but won't get to it until monday. It is definitely an elevatr thing. I am trying to switch away from proj4 strings to wkt. I must have missed something with the raster input locations. Stay tuned!

On Fri, Jan 8, 2021 at 8:31 PM Colin J. Carlson notifications@github.com wrote:

Ah I spoke too soon I'm sorry to say. Here's a reproducible example:

library(elevatr) library(raster)

sessionInfo()

rodents <- raster('RodentRichness.tif') rodents <- trim(rodents) elev <- elevatr::get_elev_raster(rodents,z=6)

This works on my laptop but not the RStudio server I'm working on, which has switched to a new error:

elev <- elevatr::get_elev_raster(rodents,z=6) Error in loc_check(locations, prj) : Please supply a valid WKT string. In addition: Warning messages: 1: In sp::wkt(locations) : no wkt comment 2: In sp::wkt(locations) : no wkt comment 3: In sp::wkt(locations) : no wkt comment

sessionInfo() on server:

sessionInfo() R version 4.0.2 (2020-06-22) Platform: x86_64-pc-linux-gnu (64-bit) Running under: Ubuntu 16.04.7 LTS

Matrix products: default BLAS: /usr/lib/libblas/libblas.so.3.6.0 LAPACK: /usr/lib/lapack/liblapack.so.3.6.0

locale: [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C LC_TIME=en_US.UTF-8 [4] LC_COLLATE=en_US.UTF-8 LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8 [7] LC_PAPER=en_US.UTF-8 LC_NAME=C LC_ADDRESS=C [10] LC_TELEPHONE=C LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

attached base packages: [1] stats graphics grDevices utils datasets methods base

other attached packages: [1] raster_3.4-5 sp_1.4-4 elevatr_0.3.3

loaded via a namespace (and not attached): [1] compiler_4.0.2 rgdal_1.5-19 tools_4.0.2 Rcpp_1.0.5 codetools_0.2-18 grid_4.0.2 [7] lattice_0.20-41

sessionInfo() on laptop:

sessionInfo() R version 4.0.3 (2020-10-10) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 10 x64 (build 18363)

Matrix products: default

locale: [1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United States.1252 [3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C [5] LC_TIME=English_United States.1252

attached base packages: [1] stats graphics grDevices utils datasets methods base

other attached packages: [1] raster_3.4-5 sp_1.4-4 elevatr_0.3.3

loaded via a namespace (and not attached): [1] compiler_4.0.3 tools_4.0.3 Rcpp_1.0.5 codetools_0.2-16 [5] grid_4.0.3 lattice_0.20-41

rodent file to play with:

Layers.zip https://github.com/jhollist/elevatr/files/5789957/Layers.zip

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/jhollist/elevatr/issues/38#issuecomment-757072852, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABJPYSY7ZFHSXUZQA2GYHWDSY6WXJANCNFSM4VPVZUHQ .

-- Jeffrey W. Hollister email: jeff.w.hollister@gmail.com cell: 401 556 4087 https://jwhollister.com

[cjcarlson]

Hey Jeff! Just pinging this, not with any urgency (because we do currently all live in hell) but just as a curiosity :)

[jhollist]

New version on CRAN should fix this. I was able to run your example without any issues. Let me know if it works for you too. I think I have a better handle on all the coordinate reference system changes. It should work for most (I won't say all) CRS use cases.

[cjcarlson]

Hey man - I think I found the issue. If I include a prj argument, I can get to a much later error message:

> elev <- elevatr::get_elev_raster(rodents, prj = crs(rodents), z=6)
Error in sp::CRS(SRS_string = prj$wkt) : 
  no arguments in initialization list
In addition: Warning messages:
1: In sp::wkt(locations) : no wkt comment
2: In sp::wkt(locations) : no wkt comment
3: In sp::wkt(locations) : no wkt comment
4: In sp::wkt(locations) : no wkt comment
5: In sp::wkt(locations) : no wkt comment
6: In sp::wkt(locations) : no wkt comment

So I dug through the code, and got to the offending step:

> prj$wkt
[1] "GEOGCS[\"WGS 84\",\n    DATUM[\"WGS_1984\",\n        SPHEROID[\"WGS 84\",6378137,298.257223563,\n            AUTHORITY[\"EPSG\",\"7030\"]],\n        AUTHORITY[\"EPSG\",\"6326\"]],\n    PRIMEM[\"Greenwich\",0,\n        AUTHORITY[\"EPSG\",\"8901\"]],\n    UNIT[\"degree\",0.0174532925199433,\n        AUTHORITY[\"EPSG\",\"9122\"]],\n    AUTHORITY[\"EPSG\",\"4326\"]]"
> sp::CRS(SRS_string = prj$wkt)
Error in sp::CRS(SRS_string = prj$wkt) : 
  no arguments in initialization list

I think this happens because my server is running GDAL 2.2.2 but the documentation says "only used when rgdal is built with PROJ >= 6 and GDAL >= 3; a valid WKT string or SRS definition such as "EPSG:4326" or "ESRI:102761"" - I think we're using old GDAL, it doesn't pick up the SRS_string argument, it thinks the function is empty, bam. I'm gonna ask sysadmins to update and will get back to you but in the meantime I think you might be able to sleep soundly

[cjcarlson]

they're also running PROJ runtime: Rel. 4.9.2, 08 September 2015 so... yeah. bam.

[jhollist]

That'd do it!

[jhollist]

@cjcarlson did they update PROJ on your server? Just curious if that fixed it and if I can close the issue.

[cjcarlson]

nope, it's going to take a thousand years apparently! It's definitely that though. Thanks for all your help!