ncss-tech / soilDB

soilDB: Simplified Access to National Cooperative Soil Survey Databases
http://ncss-tech.github.io/soilDB/
79 stars 19 forks source link

alignment of WCS results with authoritative grid topology #313

Closed brownag closed 8 months ago

brownag commented 8 months ago

As a proof of concept I use terra::align() to adjust results based on known dimensions and extent of the complete CONUS/EPSG:5070 30m data sources.

brownag commented 8 months ago

It is awkward that HI and PR have a non-integer origin. If that is correct @dylanbeaudette I think this should resolve #307 for the time being

dylanbeaudette commented 8 months ago

Meant to post this here:

OK, new HI and PR grids created, and pushed to the WCS servers. Here are the source grid specs:

HI

class       : SpatRaster 
dimensions  : 12441, 17193, 1  (nrow, ncol, nlyr)
resolution  : 30, 30  (x, y)
extent      : 56992, 572782, 8585, 381815  (xmin, xmax, ymin, ymax)
coord. ref. : NAD83(PA11) / Hawaii zone 1 (EPSG:6628) 
source      : hi-mukey-grid.tif 
name        : mukey_int 
min value   :    467795 
max value   :   3172064 

PR

class       : SpatRaster 
dimensions  : 2229, 9608, 1  (nrow, ncol, nlyr)
resolution  : 30, 30  (x, y)
extent      : 39905, 328145, 208815, 275685  (xmin, xmax, ymin, ymax)
coord. ref. : NAD83 / Puerto Rico & Virgin Is. (EPSG:32161) 
source      : pr-mukey-grid.tif 
name        : mukey_int 
min value   :    326380 
max value   :   3349465
brownag commented 8 months ago

Great! I now get:

library(soilDB )
library(sf)
#> Linking to GEOS 3.11.2, GDAL 3.6.2, PROJ 9.2.0; sf_use_s2() is TRUE
library(terra)
#> terra 1.7.55

# paste your BBOX text here
bb <- '-65.7741 18.1711,-65.7741 18.3143,-65.5228 18.3143,-65.5228 18.1711,-65.7741 18.1711'

# convert WKT string -> sfc geometry
wkt <- sprintf('POLYGON((%s))', bb)
x <- st_as_sfc(wkt, crs = 4326)

# query WCS
mu <- mukey.wcs(x, db = 'pr_ssurgo')

# make missing data (NA; the ocean) blue
plot(
  mu,
  legend = FALSE,
  axes = FALSE,
  main = attr(mu, 'layer name'),
  colNA = 'royalblue'
)

mu
#> class       : SpatRaster 
#> dimensions  : 532, 888, 1  (nrow, ncol, nlyr)
#> resolution  : 30, 30  (x, y)
#> extent      : 269705, 296345, 237495, 253455  (xmin, xmax, ymin, ymax)
#> coord. ref. : NAD83 / Puerto Rico & Virgin Is. (EPSG:32161) 
#> source(s)   : memory
#> varname     : file40245a7534fd 
#> categories  : mukey 
#> name        :   mukey 
#> min value   :  326853 
#> max value   : 3349461

# paste your BBOX text here
bb <- '-159.7426 21.9059,-159.7426 22.0457,-159.4913 22.0457,-159.4913 21.9059,-159.7426 21.9059'

# convert WKT string -> sfc geometry
wkt <- sprintf('POLYGON((%s))', bb)
x <- st_as_sfc(wkt, crs = 4326)

# query WCS
mu <- mukey.wcs(x, db = 'hi_ssurgo')

# make NA (the ocean) blue
plot(
  mu,
  legend = FALSE,
  axes = FALSE,
  main = attr(mu, 'layer name'),
  colNA = 'royalblue'
)

mu
#> class       : SpatRaster 
#> dimensions  : 540, 881, 1  (nrow, ncol, nlyr)
#> resolution  : 30, 30  (x, y)
#> extent      : 61342, 87772, 345515, 361715  (xmin, xmax, ymin, ymax)
#> coord. ref. : NAD83(PA11) / Hawaii zone 1 (EPSG:6628) 
#> source(s)   : memory
#> varname     : file402415803d8 
#> categories  : mukey 
#> name        :  mukey 
#> min value   : 467795 
#> max value   : 467929