Add interface to various RAP APIs

[brownag]

An extension of #10 and #11 to build out reusable machinery to cover tabular API endpoints other than 16-day production:

• coverV3(feat: add get_rap_table() unified interface to tabular data API endpoints #20)

• coverMeteorologyV3 (feat: add get_rap_table() unified interface to tabular data API endpoints #20)

• productionV3 (feat: add get_rap_table() unified interface to tabular data API endpoints #20)

• Add README section on get_rap_table()

• Add vignette section on get_rap_table()

Also, create a helper function for constructing URLs for map tile endpoints

• Map Tiles
  • usda-rap-tiles-cover-v3/masked/{vegetation-type}/{year}/{z}/{x}/{y}.png
  • usda-rap-tiles-biomass-v3/masked/{vegetation-type}/{year}/{z}/{x}/{y}.png
  • usda-rap-tiles-biomass-v3/unmasked/{vegetation-type}/{year}/{z}/{x}/{y}.png

https://rangelands.app/support/71-api-documentation

[brownag]

The function naming style for these tabular API endpoint functions will be get_rap_*_table() where * reflects the endpoint but not the version. E.g. get_rap_production16day_table(). This makes the naming unambiguous where we support both raster and tabular versions of the same data product. Any future versions of the API will be handled as additional argument version to the function (similar to how it is currently handled in get_rap())

[christinarwinters]

Hi there! I am struggling with the get_rap_production16day_table() continually failing and I can't figure out why. For reference, I am using gps locations to get 16 day tables for the years 2017-2024. The code runs smoothly for 2017-2019 and then starts failing at 2020 and onward. I am running it by using:

test7 <- get_rap_production16day_table(

• covs_2023,
• year = 2023,
• mask = FALSE,
• nodata_flag = NA_real_
• )

I initially had it in a loop to run through multiple years but then separated it by years to try to isolate the failure point.

The error I keep encountering is as follows:

Error in $<-.data.frame(*tmp*, "feature", value = 689L) :
replacement has 1 row, data has 0

I'm sure I'm just missing something but figured I'd ask here in case there is any insight.

Thanks!
Christina

[brownag]

Hi @christinarwinters, thanks for letting me know of the trouble.

While I can't say for sure without seeing the specific data you are using, I think that one or more of the coordinates in your input is an EMPTY geometry that returns no results (empty data matrix after removal of headers). The 689th point in covs_2023 would be one to look at.

I think I can minimally reproduce the error you are having as follows:

library(rapr)
aoi <- terra::vect("POINT EMPTY", crs = "EPSG:4326")
res <- get_rap_production16day_table(aoi, year = 2023, mask = FALSE)
#> Loading required namespace: jsonlite
#> Request failed [500]. Retrying in 1.9 seconds...
#> Request failed [500]. Retrying in 1.7 seconds...
#> Error in `$<-.data.frame`(`*tmp*`, "feature", value = 1L): replacement has 1 row, data has 0

You can use terra is.na() on AOI to get a vector where TRUE values are points with NA coordinates (empty geometry)
Quickest solution for you is to implement that filtering for yourself prior to passing covs_2023 to rapr. I can also make the code robust to this possible error so it wont mess up multi-point, multi-year runs.

In all of my testing invalid coordinates/coordinates out of the domain return -99 (NoData) for the requested interval, but I did not handle the possibility of empty geometries, which it appears the API will return error 500 (internal server error) from.

I will open up a separate issue for this, but it is relevant to this issue as well as we will need to have the same handling in other API endpoints. Thanks again for bringing it up

[christinarwinters]

Hi there!

Thank you so much for getting back to me! That was my first thought as well but it changes the failure point each time I try it.

test8 <- get_rap_production16day_table(

• covs_2024,
• year = 2024,
• mask = FALSE,
• nodata_flag = NA_real_
• )
  Error in $<-.data.frame(*tmp*, "feature", value = 21L) :
  replacement has 1 row, data has 0

test8 <- get_rap_production16day_table(

• covs_2024,
• year = 2024,
• mask = FALSE,
• nodata_flag = NA_real_
• )
  Error in $<-.data.frame(*tmp*, "feature", value = 17L) :
  replacement has 1 row, data has 0

test8 <- get_rap_production16day_table(

• covs_2024,
• year = 2024,
• mask = FALSE,
• nodata_flag = NA_real_
• )
  Request failed [500]. Retrying in 1.4 seconds...
  Error in $<-.data.frame(*tmp*, "feature", value = 27625L) :
  replacement has 1 row, data has 0

And I ran is.na() as well to check in case I was just not seeing empties but there isn't any. It's a rather big request with ~45,000 points in some of the AOI requests. I'm kind of stumped on what's going on at this point. I greatly appreciate you helping me work through this!

Thank you so much!

[brownag]

Ah, I see. I figured you might be running a lot of points. Point query APIs like this are often not ideal for high throughput applications.

I just merged #16 which adds some handling to get_rap_production16day_table() in the event that after retry we continue to get an empty response. The changes in there should fix empty results arising from error 500 in general--even if not due to empty geometry. That feature will be omitted from the result and will not throw a hard error.

Please install latest rapr from GitHub, and when you run you will see a warning with any feature ID that returned no data.

---

Assuming the coordinates are valid, the core issue seems to be server-side: under some conditions server does not generate a valid JSON response. This could be due to malformed input (as was the case with empty geoms), or some other internal failure.

There may be rate-limiting or resource limit you are hitting with large number of requests. #14 is a similar case where we run into issues with network intermittency in getting remote files with GDAL. This issue can be difficult to debug as it depends on network, pattern of use, user/project/server load, etc. and is not a "defined" behavior for the API.

If after #16 you are getting many feature IDs throwing warnings/NULL results, try to batch your AOIs into smaller chunks (say, ~1000 points each), with a brief period of waiting between each chunk, and then combine the chunked results. If doing that resolves your issue, then that might suggest there are implicit rate limitations we could better address in the function when user input is large.

I have an issue open (#13) to implement get_rap()-like (gridded/raster) functionality with 16-day data. This would require using an interface to the Python Earth Engine API (i.e. rgeedim, rgee) or GDAL to download the data. I have not benchmarked this against spatially exhaustive point queries, but my guess is it would be faster to work directly with the rasters.