Geomorphology doesn't seem to work for gridded data

[jaseeverett]

Running the Bermuda example in the help fails

library(oceandatr)
# Grab EEZ data first 
bermuda_eez <- get_area(area_name = "Bermuda", mregions_column = "territory1")
# Get geomorphological features in spatial_grid
planning_grid <- get_grid(area_polygon = bermuda_eez, projection_crs = '+proj=laea +lon_0=-64.8108333 +lat_0=32.3571917 +datum=WGS84 +units=m +no_defs', resolution = 5000)
geomorph_gridded <- get_geomorphology(spatial_grid = planning_grid)
#> Error in sf_to_grid(dat, spatial_grid, matching_crs, name, feature_names, : No data in grid.

^{Created on 2024-04-02 with reprex v2.1.0}

[jflowernet]

Thanks Jase. I don't get an error. I think it might be due to recent updates to spatialgridr and oceandatr: can you try installing the latest versions of both packages and try again?

[jaseeverett]

No that doesn't seem to be the problem on my end. Still getting the same error.

But note, it doesn't look like you increasing the version of spatialgridr when you submit changes. If spatialgridr and oceandatr are dependent on certain versions, it would be good to set minimum versions in the DESCRIPTION and also increment the version number each time you push. I had to do a force=TRUE to get spatialdridr to install because it thinks it's the same version.
Also I am getting a warning about some of your data products needing at least R > 3.5, so it would be good to set that dependency in the DESCRIPTION as well.

I'll try and do some more debugging tomorrow. Gotta prepare for a meeting now

[jflowernet]

I re-installed both packages from Github and still don't get any errors. Any extra info you can give would be great.

Noted about the dependencies and versioning. I'll add the R>=3.5 dependency to both. Using devtools::install_github() checks the package is up to date with the latest commit, e.g. I get this if I have exactly the same version install:
Skipping install of 'spatialgridr' from a github remote, the SHA1 (8991bac3) has not changed since last install. Use force = TRUE to force installation

So you shouldn't have to use force = TRUE

[jaseeverett]

Fair call about the force = TRUE.
But also add a dependency to spatialgridr version >= in oceandatr. Otherwise R won't force an update to spatialdatr when oceandatr is updated.

[jaseeverett]

Both out2 and out3 fail for me. What about you? Perhaps its because there is no data and there is no fall back?

area <- oceandatr::get_area("Maldives")
PUs <- spatialgridr::get_grid(area_polygon = area, projection_crs = "ESRI:54009")
dat <- oceandatr::get_knolls(area_polygon = area)
out <- spatialgridr::get_data_in_grid(spatial_grid = PUs, dat = dat, apply_cutoff = FALSE)

# Setup area with completely different limits
area2 <- oceandatr::get_area("Brazil")
PUs2 <- spatialgridr::get_grid(area_polygon = area2, projection_crs = "ESRI:54009")

# Try and use area2 to cross with data from dat
out2 <- spatialgridr::get_data_in_grid(spatial_grid = PUs2, dat = dat, apply_cutoff = TRUE, cutoff = 0.001)
out3 <- spatialgridr::get_data_in_grid(spatial_grid = PUs2, dat = dat, apply_cutoff = FALSE)

[jflowernet]

Yep definitely fail for out2 and 3, but that's because it is trying to intersect data from the Maldives with the Brazil PUs! You should get a slightly helpful error message (this is what I get):

Error in sf_to_grid(dat, spatial_grid, matching_crs, name, feature_names,  : 
  No data in grid.
In addition: Warning message:
attribute variables are assumed to be spatially constant throughout all geometries 

That being said, when I try

dat2 <- oceandatr::get_knolls(area2)
out4 <- spatialgridr::get_data_in_grid(spatial_grid = PUs2, dat = dat2, apply_cutoff = TRUE, cutoff = 0.001)

I get an error because unique(sf::st_geometry_type(dat2)) returns POLYGON MULTIPOLYGON. I guess I could add a catch to get_data_in_grid to check for this and st_cast(to = "MULTIPOLYGON") if there are both polygon types present? But I think it's not totally unreasonable to expect the input to be a single geometry type!

[jaseeverett]

I'll test the 2nd error you highlight. sf should be able to deal with that.

But the first one (my example I sent) shouldn't error I don't think. It should return the PUs with 0's in each row. Imagine you have global seamount data with a small planning area and there are no seamounts within the planning area. It shouldn't error but return 0's.

Perhaps I should be creating that exact scenario, but I didn't have a chance to. Basically if there is no coverage it should return no coverage, not an error. Does that make sense?

[jflowernet]

Yeah, so that's an interesting point. I actually changed the package to return error when there is no data in the area (area polygon or spatial grid). My thinking was that the user should be made aware of this: if they just get an empty grid/ polygon, surely you'd be thinking "where's my data"?
However, this has caused me some problems, and I suspect the same for you, when you are trying to loop through data using get_data_in_grid. I had to add tryCatch to some of the oceandatr functions to avoid it stopping as soon as the first no data grid is found.

I'd be interested to get your thoughts though.

[jaseeverett]

I think you should return a warning or display a message, but then return the PUs with 0s. It stuffs up my pipes or scripts if it errors. For example in shinyplanr, I throw all the data at it, but it doesn't mean there will be data in all PUs. It's the same if the data doesn't meet the cutoff. It returns a 0, not an error.

How would it work with your get_features() script didn't have a seamount, but had a pelagic regionalisation. It would error I assume, despite the fact some data exists.

Happy to chat in more detail after the meeting tomorrow if you want.