search

mapme-initiative / mapme.biodiversity

Efficient analysis of spatial biodiversity datasets for global portfolios
https://mapme-initiative.github.io/mapme.biodiversity/dev
GNU General Public License v3.0
33 stars 7 forks source link

CRS issue with frits_et_al #300

Closed fBedecarrats closed 3 months ago

fBedecarrats commented 3 months ago

The following generates warnings including error messages:

library(tidyverse)
library(sf)
library(mapme.biodiversity) # This is CRAN version
library(geodata)

mada <- gadm("Madagascar", level = 0, resolution = 2, path = ".")
mada_dt <- mada %>%
  st_as_sf() %>%
  st_cast("POLYGON")

st_crs(mada)
mada_dt <- mada_dt %>%
  get_resources(get_fritz_et_al()) %>%
  calc_indicators(calc_deforestation_drivers())

Here is the output:

R version 4.4.0 (2024-04-24) -- "Puppy Cup"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

  Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

Connected to your session in progress, last started 2024-Jul-05 08:47:47 UTC (2 hours ago)
> install.packages("mapme.biodiversity")
Installing package into ‘/usr/local/lib/R/site-library’
(as ‘lib’ is unspecified)
also installing the dependencies ‘listenv’, ‘parallelly’, ‘future’, ‘globals’, ‘furrr’

trying URL 'https://packagemanager.posit.co/cran/__linux__/jammy/latest/src/contrib/listenv_0.9.1.tar.gz'
Content type 'binary/octet-stream' length 106238 bytes (103 KB)
==================================================
downloaded 103 KB

trying URL 'https://packagemanager.posit.co/cran/__linux__/jammy/latest/src/contrib/parallelly_1.37.1.tar.gz'
Content type 'binary/octet-stream' length 361300 bytes (352 KB)
==================================================
downloaded 352 KB

trying URL 'https://packagemanager.posit.co/cran/__linux__/jammy/latest/src/contrib/future_1.33.2.tar.gz'
Content type 'binary/octet-stream' length 656257 bytes (640 KB)
==================================================
downloaded 640 KB

trying URL 'https://packagemanager.posit.co/cran/__linux__/jammy/latest/src/contrib/globals_0.16.3.tar.gz'
Content type 'binary/octet-stream' length 108405 bytes (105 KB)
==================================================
downloaded 105 KB

trying URL 'https://packagemanager.posit.co/cran/__linux__/jammy/latest/src/contrib/furrr_0.3.1.tar.gz'
Content type 'binary/octet-stream' length 1019440 bytes (995 KB)
==================================================
downloaded 995 KB

trying URL 'https://packagemanager.posit.co/cran/__linux__/jammy/latest/src/contrib/mapme.biodiversity_0.8.0.tar.gz'
Content type 'binary/octet-stream' length 782384 bytes (764 KB)
==================================================
downloaded 764 KB

* installing *binary* package ‘listenv’ ...
* DONE (listenv)
* installing *binary* package ‘parallelly’ ...
* DONE (parallelly)
* installing *binary* package ‘globals’ ...
* DONE (globals)
* installing *binary* package ‘future’ ...
* DONE (future)
* installing *binary* package ‘furrr’ ...
* DONE (furrr)
* installing *binary* package ‘mapme.biodiversity’ ...
* DONE (mapme.biodiversity)

The downloaded source packages are in
    ‘/tmp/RtmpjK29Zs/downloaded_packages’
> mada <- gadm("Madagascar", level = 0, resolution = 2)
Error in gadm("Madagascar", level = 0, resolution = 2) : 
  could not find function "gadm"
> library(tidyverse)
── Attaching core tidyverse packages ──────────────────────────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.4     ✔ readr     2.1.5
✔ forcats   1.0.0     ✔ stringr   1.5.1
✔ ggplot2   3.5.1     ✔ tibble    3.2.1
✔ lubridate 1.9.3     ✔ tidyr     1.3.1
✔ purrr     1.0.2     
── Conflicts ────────────────────────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package to force all conflicts to become errors
> library(sf)
Linking to GEOS 3.10.2, GDAL 3.4.1, PROJ 8.2.1; sf_use_s2() is TRUE
> library(mapme)
Error in library(mapme) : there is no package called ‘mapme’
> mada_dt < mada %>%
+   st_cast("POLYGON")
Error: object 'mada_dt' not found
> mada_dt <- mada %>%
+   st_cast("POLYGON")
Error: object 'mada' not found
> library(mapme.biodiversity)
> library(geodata)
Loading required package: terra
terra 1.7.78

Attaching package: ‘terra’

The following object is masked from ‘package:tidyr’:

    extract

> mada <- gadm("Madagascar", level = 0, resolution = 2)
Error: path is missing
> mada <- gadm("Madagascar", level = 0, resolution = 2, path = ".")
trying URL 'https://geodata.ucdavis.edu/gadm/gadm4.1/pck/gadm41_MDG_0_pk_low.rds'
Content type 'unknown' length 84541 bytes (82 KB)
==================================================
downloaded 82 KB
> st_crs(mada)
Coordinate Reference System:
  User input: WGS 84 
  wkt:
GEOGCRS["WGS 84",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["geodetic latitude (Lat)",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["geodetic longitude (Lon)",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    ID["EPSG",4326]]
> mada_dt <- mada %>%
+   st_to_sf() %>%
+   st_cast("POLYGON")
Error in st_to_sf(.) : could not find function "st_to_sf"
> mada_dt <- mada %>%
+   st_as_sf() %>%
+   st_cast("POLYGON")
Warning message:
In st_cast.sf(., "POLYGON") :
  repeating attributes for all sub-geometries for which they may not be constant
> mada_dt <- mada_dt %>%
+   get_resources(get_fritz_et_al()) %>%
+   calc_indicators(calc_deforestation_drivers())
Found a column named 'assetid'. Overwritting its values with a unique identifier.
There were 50 or more warnings (use warnings() to see the first 50)
> warnings()
Warning messages:
1: In check_namespace("progressr", error = FALSE) :
  R package 'progressr' required.
Please install via `install.packages('progressr')`FALSE
2: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
3: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
4: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
5: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
6: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
7: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
8: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
9: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
10: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
11: In CPL_gdalinfo(if (missing(source)) character(0) else source,  ... :
  GDAL Error 1: Point outside of projection domain
12: In check_namespace("progressr", error = FALSE) :
  R package 'progressr' required.
Please install via `install.packages('progressr')`FALSE
13: [mask] CRS do not match
14: [mask] CRS do not match
15: [mask] CRS do not match
16: [mask] CRS do not match
17: [mask] CRS do not match
18: [mask] CRS do not match
19: [mask] CRS do not match
20: [mask] CRS do not match
21: [mask] CRS do not match
22: [mask] CRS do not match
23: [mask] CRS do not match
24: [mask] CRS do not match
25: [mask] CRS do not match
26: [mask] CRS do not match
27: [mask] CRS do not match
28: [mask] CRS do not match
29: [mask] CRS do not match
30: [mask] CRS do not match
31: [mask] CRS do not match
32: [mask] CRS do not match
33: [mask] CRS do not match
34: [mask] CRS do not match
35: [mask] CRS do not match
36: [mask] CRS do not match
37: [mask] CRS do not match
38: [mask] CRS do not match
39: [mask] CRS do not match
40: [mask] CRS do not match
41: [mask] CRS do not match
42: [mask] CRS do not match
43: [mask] CRS do not match
44: [mask] CRS do not match
45: [mask] CRS do not match
46: [mask] CRS do not match
47: [mask] CRS do not match
48: [mask] CRS do not match
49: [mask] CRS do not match
50: [mask] CRS do not match

And the results are all in 0: 

> > mada_dt %>%
+   group_by(COUNTRY, datetime, variable) %>%
+   st_drop_geometry() %>%
+   summarise(n = n(),
+             total = sum(value, na.rm = FALSE))
`summarise()` has grouped output by 'COUNTRY', 'datetime'. You can override using the `.groups`
argument.
# A tibble: 10 × 5
# Groups:   COUNTRY, datetime [1]
   COUNTRY    datetime            variable                          n total
   <chr>      <dttm>              <chr>                         <int> <dbl>
 1 Madagascar 2008-01-01 00:00:00 commercial_agriculture          416     0
 2 Madagascar 2008-01-01 00:00:00 commercial_oil_palm             416     0
 3 Madagascar 2008-01-01 00:00:00 managed_forests                 416     0
 4 Madagascar 2008-01-01 00:00:00 mining                          416     0
 5 Madagascar 2008-01-01 00:00:00 natural_disturbances            416     0
 6 Madagascar 2008-01-01 00:00:00 other_subsistance_agriculture   416     0
 7 Madagascar 2008-01-01 00:00:00 pasture                         416     0
 8 Madagascar 2008-01-01 00:00:00 roads                           416     0
 9 Madagascar 2008-01-01 00:00:00 shifting_cultivation            416     0
10 Madagascar 2008-01-01 00:00:00 wildfire                        416     0
goergen95 commented 3 months ago

Thanks, for catching this. We do not transform the asset to the CRS of the resource in calc_deforestation_area(). I'll push a fix soon. In the meantime, this is how you can get to the raster:

library(sf)
library(terra)
library(mapme.biodiversity) # This is CRAN version
library(geodata)

mada <- st_as_sf(gadm("Madagascar", level = 0, resolution = 2, path = "."))
get_resources(mada, get_fritz_et_al())
fritz_et_al <- prep_resources(mada, resources = "fritz_et_al")[[1]]
fBedecarrats commented 3 months ago

That was quick! Thank you so much!