st_as_stars method for ncdfgeom object Fixes #65

DESCRIPTION

@@ -1,15 +1,15 @@
 Package: ncdfgeom
 Type: Package
 Title: 'NetCDF' Geometry and Time Series
-Version: 1.0.0
+Version: 1.1.0
 Date: 2019-06-05
 Authors@R: c(person("David", "Blodgett", role = c("aut", "cre"), 
                     email = "dblodgett@usgs.gov"),
              person("Luke", "Winslow", role = "ctb"))
 Description: Tools to create time series and geometry 'NetCDF' files.
 URL: https://code.usgs.gov/water/ncdfgeom
 BugReports: https://github.com/USGS-R/ncdfgeom/issues
-Imports: RNetCDF, ncmeta, sf, dplyr, methods
+Imports: RNetCDF, ncmeta, sf, dplyr, methods, stars
 Depends:
   R (>= 3.0)
 Suggests: testthat, knitr, rmarkdown, pkgdown, tidyverse, sp, geoknife, ncdf4, jsonlite

NAMESPACE

@@ -1,5 +1,6 @@
 # Generated by roxygen2: do not edit by hand
 
+S3method(st_as_stars,ncdfgeom)
 export(read_attribute_data)
 export(read_geometry)
 export(read_timeseries_dsg)
@@ -38,4 +39,5 @@ importFrom(sf,st_set_geometry)
 importFrom(sf,st_sf)
 importFrom(sf,st_sfc)
 importFrom(sf,st_zm)
+importFrom(stars,st_as_stars)
 importFrom(stats,setNames)

R/st_as_stars.R

@@ -0,0 +1,43 @@
+#' Convert ncdfgeom object into stars object.
+#' @importFrom stars st_as_stars
+#' @param .x Object of class ncdfgeom as returned by read_timeseries_dsg.
+#' @param ... not used.
+#' @param sf_geometry sf data.frame with geometry and attributes to be added to stars object. 
+#' Must have same number of rows as timeseries instances.
+#' @name st_as_stars
+#' @export
+#' 
+st_as_stars.ncdfgeom <- function(.x, ..., sf_geometry = NA) {
+  crs <- st_crs(4326)$proj4string
+  ts_points <- data.frame(X = .x$lons, Y = .x$lats, Z = .x$alts)
+  ts_points <- sf::st_as_sf(ts_points, coords = c("X", "Y", "Z"), crs = crs)
+
+  data <- .x$data_frames[[1]]
+  # data[["T"]] <- .x$time
+
+  gdim <- stars:::create_dimension(from = 1, to = length(.x$lats), 
+                                   refsys = crs, point = TRUE, 
+                                   values = ts_points$geometry)
+  tdim <- stars:::create_dimension(from = 1, to = length(.x$time), 
+                                   refsys = "POSIXct", point = FALSE, 
+                                   values = as.POSIXct(.x$time))
+  dim <- list(time = tdim, points = gdim)
+
+  if("sf" %in% class(sf_geometry)) {
+    if(length(gdim$values) != length(st_geometry(sf_geometry))) 
+      stop("geometry must be same length as instance dimension of timeseries")
+
+    is_point <- any(grepl("point", class(st_geometry(sf_geometry)), ignore.case = TRUE))
+
+    sf_dim <- stars:::create_dimension(from = 1, to = length(gdim$values),
+                                       refsys = st_crs(sf_geometry)$proj4string, 
+                                       point = is_point, is_raster = FALSE,
+                                       values = st_geometry(sf_geometry))
+
+    dim <- c(dim, list(geometry = sf_dim))
+  }
+
+  stars:::st_stars(x = setNames(list(as.matrix(.x$data_frames[[1]])), 
+                                .x$varmeta[[1]]$name), 
+                   dimensions =  stars:::create_dimensions(dim))
+}

tests/testthat/helper-functions.R

@@ -97,5 +97,36 @@ get_sample_timeseries_data <- function() {
               lons = lons,
               lats = lats,
               alts = alts,
-              units = units))
+              units = units,
+              geom = yahara))
+}
+
+get_test_ncdf_object <- function(nc_file = tempfile()) {
+  nc_summary<-'test summary'
+  nc_date_create<-'2099-01-01'
+  nc_creator_name='test creator'
+  nc_creator_email='test@test.com'
+  nc_project='testthat ncdfgeom'
+  nc_proc_level='just a test no processing'
+  nc_title<-'test title'
+  global_attributes<-list(title = nc_title, summary = nc_summary, date_created=nc_date_create, 
+                          creator_name=nc_creator_name,creator_email=nc_creator_email,
+                          project=nc_project, processing_level=nc_proc_level)
+
+  test_data <- get_sample_timeseries_data()
+
+  testnc<-write_timeseries_dsg(nc_file, 
+                               names(test_data$var_data), 
+                               test_data$lats, test_data$lons, 
+                               as.character(test_data$time), 
+                               test_data$var, 
+                               test_data$alts, 
+                               data_unit=test_data$units,	
+                               data_prec='double',
+                               data_metadata=test_data$meta,
+                               attributes=global_attributes)
+
+  test_nc <- write_geometry(nc_file, test_data$geom, variables = test_data$meta$name)
+
+  list(ncdfgeom = read_timeseries_dsg(nc_file), sf = read_geometry(nc_file))
 }

tests/testthat/test_read-write_timeseries_dsg.R

@@ -2,6 +2,7 @@ context("orthogonal netcdf timeseries")
 
 test_that("Create basic DSG file", {
 
+  # NOTE: this code has been moved to helper files but was left here to not mess with it.
   nc_file<-tempfile()
   nc_summary<-'test summary'
   nc_date_create<-'2099-01-01'

tests/testthat/test_st_as_stars.R

@@ -0,0 +1,27 @@
+context("st_as_stars tests")
+
+test_that("basic st_as_stars", {
+  test_list <- get_test_ncdf_object()
+
+  stars_obj <- st_as_stars(test_list$ncdfgeom)
+
+  expect_s3_class(stars_obj, "stars")
+
+  dim <- stars::st_dimensions(stars_obj)
+  expect_equal(sf::st_crs(dim$points$refsys), sf::st_crs(4326))
+  expect_equal(dim$time$refsys, "POSIXct")
+
+  expect_s3_class(dim$points$values, "sfc_POINT")
+
+  expect_true(dim$points$point)  
+
+  stars_obj <- st_as_stars(test_list$ncdfgeom, sf_geometry = test_list$sf)
+
+  dim <- stars::st_dimensions(stars_obj)
+  expect_equal(sf::st_crs(dim$geometry$refsys), sf::st_crs(test_list$sf))
+
+  expect_s3_class(dim$geometry$values, "sfc_POLYGON")
+
+  expect_false(dim$geometry$point)  
+
+})