r-spatial / stars

Spatiotemporal Arrays, Raster and Vector Data Cubes
https://r-spatial.github.io/stars/
Apache License 2.0
558 stars 94 forks source link

A better (efficient) way to convert a list to stars object with time dimension having only one value #691

Closed alexyshr closed 2 months ago

alexyshr commented 3 months ago

Dear All,

I need to convert the list l1 (which contains one element) into a stars object with one attribute (any name) and three dimensions (x, y, and time). The time dimension must have only one value equal to t1. The result should be equivalent to the object r1 (as shown at the end of this code). I have already solved the issue in a way that seems to overload the process (I can leave it that way if there is no better way to do it).

As this process is done inside a loop, I need the r1 object always to have a third dimension called time. Note that with more than one element in the list, there is no issue with creating the time dimension (as shown with r2)

A possible bug?: Why the name of the last dimension of r3 is values instead of time. See the line r3 <- merge(r3, name = "time", values = c(t1, t1 + 1)). It seems the parameter **name** of **merge** is not working.

library(stars)
#> Loading required package: abind
#> Loading required package: sf
#> Linking to GEOS 3.12.1, GDAL 3.8.4, PROJ 9.3.1; sf_use_s2() is TRUE

#######################################################
tif <- system.file("tif/lc.tif", package = "stars")
x <- read_stars(tif)
l1 <- list(x)
t1 <- as.Date("2018-07-31")
names(l1) <- c(t1)
# I need to convert the list l1 (one element)
# to stars object with
# 1 attribure (any name)
# and 3 dimensions (x, y, and time)
# the time dimension only has one value equal to t1
# the result need to be equal to the object r1 (end of this code)

########################################################
# If the list has two elements the approach will be
l2 <- list(x, x)
names(l2) <- c(t1, t1 + 1)
r2 <- do.call("c", l2)
# st_redimension transforms attributes to dimensions (similar to merge) while keeping their names:
# the new dimension will be called new_dim
# and the values will be the names of the list l2 c(t1, t1 + 1)
r2 <- st_redimension(r2)
# renanme the new_dim to time
r2 <- st_set_dimensions(r2, which = 3, name = "time")

# another faster approach is to use merge
r2 <- do.call("c", l2)
r2 <- merge(r2, name = "time")
r2
#> stars object with 3 dimensions and 1 attribute
#> attribute(s):
#>                  2018-07-31.2018-08-01 
#>  Evergreen Forest           : 912      
#>  Herbaceuous                : 540      
#>  Open Water                 : 504      
#>  Developed, Low Intensity   : 162      
#>  Developed, Medium Intensity:  96      
#>  (Other)                    : 284      
#>  NA's                       :5230      
#> dimension(s):
#>      from to  offset delta                    refsys point
#> x       1 84 3092415  3000 Albers Conical Equal Area FALSE
#> y       1 46   59415 -3000 Albers Conical Equal Area FALSE
#> time    1  2      NA    NA                        NA    NA
#>                      values x/y
#> x                      NULL [x]
#> y                      NULL [y]
#> time 2018-07-31, 2018-08-01

####################################################
# a not clean solution for one element list:
# repeating the attribute, then subseting
# the first value of time dimension
r1 <- do.call("c", l1)
r3 <- c(r1, r1)
# the name = "time" is not working
r3 <- merge(r3, name = "time", values = c(t1, t1 + 1))
r3
#> stars object with 3 dimensions and 1 attribute
#> attribute(s):
#>                     lc.tif.lc.tif.1 
#>  Evergreen Forest           : 912   
#>  Herbaceuous                : 540   
#>  Open Water                 : 504   
#>  Developed, Low Intensity   : 162   
#>  Developed, Medium Intensity:  96   
#>  (Other)                    : 284   
#>  NA's                       :5230   
#> dimension(s):
#>        from to     offset  delta                    refsys point x/y
#> x         1 84    3092415   3000 Albers Conical Equal Area FALSE [x]
#> y         1 46      59415  -3000 Albers Conical Equal Area FALSE [y]
#> values    1  2 2018-07-31 1 days                      Date    NA
#
# as the parameter name is not working
# change the last dimension to name "time" using st_set_dimensions
r3 <- st_set_dimensions(r3, 3, name = "time")
# subset the first value of time dimension
r1 <- r3[, , , 1]
r1
#> stars object with 3 dimensions and 1 attribute
#> attribute(s):
#>                     lc.tif.lc.tif.1 
#>  Evergreen Forest           : 456   
#>  Herbaceuous                : 270   
#>  Open Water                 : 252   
#>  Developed, Low Intensity   :  81   
#>  Developed, Medium Intensity:  48   
#>  (Other)                    : 142   
#>  NA's                       :2615   
#> dimension(s):
#>      from to     offset  delta                    refsys point x/y
#> x       1 84    3092415   3000 Albers Conical Equal Area FALSE [x]
#> y       1 46      59415  -3000 Albers Conical Equal Area FALSE [y]
#> time    1  1 2018-07-31 1 days                      Date    NA

Created on 2024-06-14 with reprex v2.1.0

Session info ``` r sessioninfo::session_info() #> - Session info --------------------------------------------------------------- #> setting value #> version R version 4.4.0 (2024-04-24 ucrt) #> os Windows 10 x64 (build 19045) #> system x86_64, mingw32 #> ui RTerm #> language (EN) #> collate English_United States.utf8 #> ctype English_United States.utf8 #> tz America/New_York #> date 2024-06-14 #> pandoc 2.17.1.1 @ C:/Users/500596~1/AppData/Local/Pandoc/ (via rmarkdown) #> #> - Packages ------------------------------------------------------------------- #> package * version date (UTC) lib source #> abind * 1.4-5 2016-07-21 [1] CRAN (R 4.4.0) #> class 7.3-22 2023-05-03 [2] CRAN (R 4.4.0) #> classInt 0.4-10 2023-09-05 [2] CRAN (R 4.4.0) #> cli 3.6.2 2023-12-11 [2] CRAN (R 4.4.0) #> DBI 1.2.3 2024-06-02 [1] CRAN (R 4.4.0) #> digest 0.6.35 2024-03-11 [2] CRAN (R 4.4.0) #> e1071 1.7-14 2023-12-06 [2] CRAN (R 4.4.0) #> evaluate 0.23 2023-11-01 [2] CRAN (R 4.4.0) #> fastmap 1.1.1 2023-02-24 [2] CRAN (R 4.4.0) #> fs 1.6.4 2024-04-25 [2] CRAN (R 4.4.0) #> glue 1.7.0 2024-01-09 [2] CRAN (R 4.4.0) #> htmltools 0.5.8.1 2024-04-04 [2] CRAN (R 4.4.0) #> KernSmooth 2.23-22 2023-07-10 [2] CRAN (R 4.4.0) #> knitr 1.46 2024-04-06 [2] CRAN (R 4.4.0) #> lifecycle 1.0.4 2023-11-07 [2] CRAN (R 4.4.0) #> magrittr 2.0.3 2022-03-30 [2] CRAN (R 4.4.0) #> proxy 0.4-27 2022-06-09 [2] CRAN (R 4.4.0) #> purrr 1.0.2 2023-08-10 [2] CRAN (R 4.4.0) #> R.cache 0.16.0 2022-07-21 [2] CRAN (R 4.4.0) #> R.methodsS3 1.8.2 2022-06-13 [2] CRAN (R 4.4.0) #> R.oo 1.26.0 2024-01-24 [2] CRAN (R 4.4.0) #> R.utils 2.12.3 2023-11-18 [2] CRAN (R 4.4.0) #> Rcpp 1.0.12 2024-01-09 [2] CRAN (R 4.4.0) #> reprex 2.1.0 2024-01-11 [1] CRAN (R 4.4.0) #> rlang 1.1.4 2024-06-04 [1] CRAN (R 4.4.0) #> rmarkdown 2.26 2024-03-05 [2] CRAN (R 4.4.0) #> sessioninfo 1.2.2 2021-12-06 [2] CRAN (R 4.4.0) #> sf * 1.0-16 2024-03-24 [2] CRAN (R 4.4.0) #> stars * 0.6-5 2024-04-04 [1] CRAN (R 4.4.0) #> styler 1.10.3 2024-04-07 [2] CRAN (R 4.4.0) #> units 0.8-5 2023-11-28 [2] CRAN (R 4.4.0) #> vctrs 0.6.5 2023-12-01 [2] CRAN (R 4.4.0) #> withr 3.0.0 2024-01-16 [2] CRAN (R 4.4.0) #> xfun 0.43 2024-03-25 [2] CRAN (R 4.4.0) #> yaml 2.3.8 2023-12-11 [2] CRAN (R 4.4.0) #> #> [1] C:/Users/500596972/AppData/Local/R/win-library/4.4 #> [2] C:/Program Files/R/R-4.4.0/library #> #> ------------------------------------------------------------------------------ ```
alexyshr commented 3 months ago

A way to create an additional time dimension with one value is:

r1 <- c(x, along = list(time = t1))

In my code using a list storing stars objects, I need to do an if else section when the list has only one stars object. That particular behavior of r2 <- do.call("c", l2), then r2 <- merge(r2, name = "time") when l2 has only one element was not expected.

Am I wrong expecting that the name of the dimension will be time instead of values in r3 <- merge(r3, name = "time", values = c(t1, t1 + 1))?

alexyshr commented 3 months ago

The solution is simple. The only one that is not working is r1 <- merge(r1, time = c(t1)). It can be replaced by r1 <- c(r1, along = list(time = t1))

library(stars)
#> Loading required package: abind
#> Loading required package: sf
#> Linking to GEOS 3.12.1, GDAL 3.8.4, PROJ 9.3.1; sf_use_s2() is TRUE

#######################################################
tif <- system.file("tif/lc.tif", package = "stars")
x <- read_stars(tif)
l1 <- list(x)
t1 <- as.Date("2018-07-31")
names(l1) <- c(t1)
# Convert the list l1 (one element)
# to stars object with
# 1 attribute (time)
# and 3 dimensions (x, y, and time)
# the time dimension only has one value equal to t1
r1 <- do.call("c", l1)
# this does not work:
# r1 <- merge(r1, time = c(t1))
# but this works:
r1 <- c(r1, along = list(time = t1))

########################################################
# if the list has two elements the approach will be
l2 <- list(x, x)
names(l2) <- c(t1, t1 + 1)
r2 <- do.call("c", l2)
# st_redimension transforms attributes to dimensions (similar to merge) while keeping their names:
# the new dimension will be called new_dim
# and the values will be the names of the list l2 c(t1, t1 + 1)
r2 <- st_redimension(r2)
# renanme the new_dim to time
r2 <- st_set_dimensions(r2, which = 3, name = "time")

# another faster approach is to use merge
r2 <- do.call("c", l2)
r2 <- merge(r2, name = "time")
r2
#> stars object with 3 dimensions and 1 attribute
#> attribute(s):
#>                  2018-07-31.2018-08-01 
#>  Evergreen Forest           : 912      
#>  Herbaceuous                : 540      
#>  Open Water                 : 504      
#>  Developed, Low Intensity   : 162      
#>  Developed, Medium Intensity:  96      
#>  (Other)                    : 284      
#>  NA's                       :5230      
#> dimension(s):
#>      from to  offset delta                    refsys point
#> x       1 84 3092415  3000 Albers Conical Equal Area FALSE
#> y       1 46   59415 -3000 Albers Conical Equal Area FALSE
#> time    1  2      NA    NA                        NA    NA
#>                      values x/y
#> x                      NULL [x]
#> y                      NULL [y]
#> time 2018-07-31, 2018-08-01

####################################################
# a not clean solution for one element list:
# repeating the attribute, then subseting
# the first value of time dimension
r1 <- do.call("c", l1)
r3 <- c(r1, r1)
# the name = c("time") is ignored and the
# final name of the dimension is "the_time"
r3 <- merge(r3, name = c("time"), the_time = c(t1, t1 + 1))
r3
#> stars object with 3 dimensions and 1 attribute
#> attribute(s):
#>                     lc.tif.lc.tif.1 
#>  Evergreen Forest           : 912   
#>  Herbaceuous                : 540   
#>  Open Water                 : 504   
#>  Developed, Low Intensity   : 162   
#>  Developed, Medium Intensity:  96   
#>  (Other)                    : 284   
#>  NA's                       :5230   
#> dimension(s):
#>          from to     offset  delta                    refsys point x/y
#> x           1 84    3092415   3000 Albers Conical Equal Area FALSE [x]
#> y           1 46      59415  -3000 Albers Conical Equal Area FALSE [y]
#> the_time    1  2 2018-07-31 1 days                      Date    NA

#
# change the last dimension to name "time" using st_set_dimensions
r3 <- st_set_dimensions(r3, 3, name = "time")
# subset the first value of time dimension
r1 <- r3[, , , 1]
r1
#> stars object with 3 dimensions and 1 attribute
#> attribute(s):
#>                     lc.tif.lc.tif.1 
#>  Evergreen Forest           : 456   
#>  Herbaceuous                : 270   
#>  Open Water                 : 252   
#>  Developed, Low Intensity   :  81   
#>  Developed, Medium Intensity:  48   
#>  (Other)                    : 142   
#>  NA's                       :2615   
#> dimension(s):
#>      from to     offset  delta                    refsys point x/y
#> x       1 84    3092415   3000 Albers Conical Equal Area FALSE [x]
#> y       1 46      59415  -3000 Albers Conical Equal Area FALSE [y]
#> time    1  1 2018-07-31 1 days                      Date    NA

Created on 2024-06-29 with reprex v2.1.0

Session info ``` r sessioninfo::session_info() #> - Session info --------------------------------------------------------------- #> setting value #> version R version 4.4.0 (2024-04-24 ucrt) #> os Windows 10 x64 (build 19045) #> system x86_64, mingw32 #> ui RTerm #> language (EN) #> collate English_United States.utf8 #> ctype English_United States.utf8 #> tz America/New_York #> date 2024-06-29 #> pandoc 2.17.1.1 @ C:/Users/500596~1/AppData/Local/Pandoc/ (via rmarkdown) #> #> - Packages ------------------------------------------------------------------- #> package * version date (UTC) lib source #> abind * 1.4-5 2016-07-21 [1] CRAN (R 4.4.0) #> class 7.3-22 2023-05-03 [2] CRAN (R 4.4.0) #> classInt 0.4-10 2023-09-05 [2] CRAN (R 4.4.0) #> cli 3.6.2 2023-12-11 [2] CRAN (R 4.4.0) #> DBI 1.2.3 2024-06-02 [1] CRAN (R 4.4.0) #> digest 0.6.35 2024-03-11 [2] CRAN (R 4.4.0) #> e1071 1.7-14 2023-12-06 [2] CRAN (R 4.4.0) #> evaluate 0.23 2023-11-01 [2] CRAN (R 4.4.0) #> fastmap 1.1.1 2023-02-24 [2] CRAN (R 4.4.0) #> fs 1.6.4 2024-04-25 [2] CRAN (R 4.4.0) #> glue 1.7.0 2024-01-09 [2] CRAN (R 4.4.0) #> htmltools 0.5.8.1 2024-04-04 [2] CRAN (R 4.4.0) #> KernSmooth 2.23-22 2023-07-10 [2] CRAN (R 4.4.0) #> knitr 1.46 2024-04-06 [2] CRAN (R 4.4.0) #> lifecycle 1.0.4 2023-11-07 [2] CRAN (R 4.4.0) #> magrittr 2.0.3 2022-03-30 [2] CRAN (R 4.4.0) #> proxy 0.4-27 2022-06-09 [2] CRAN (R 4.4.0) #> purrr 1.0.2 2023-08-10 [2] CRAN (R 4.4.0) #> R.cache 0.16.0 2022-07-21 [2] CRAN (R 4.4.0) #> R.methodsS3 1.8.2 2022-06-13 [2] CRAN (R 4.4.0) #> R.oo 1.26.0 2024-01-24 [2] CRAN (R 4.4.0) #> R.utils 2.12.3 2023-11-18 [2] CRAN (R 4.4.0) #> Rcpp 1.0.12 2024-01-09 [2] CRAN (R 4.4.0) #> reprex 2.1.0 2024-01-11 [1] CRAN (R 4.4.0) #> rlang 1.1.4 2024-06-04 [1] CRAN (R 4.4.0) #> rmarkdown 2.26 2024-03-05 [2] CRAN (R 4.4.0) #> sessioninfo 1.2.2 2021-12-06 [2] CRAN (R 4.4.0) #> sf * 1.0-17 2024-06-25 [1] Github (r-spatial/sf@3705509) #> stars * 0.6-6 2024-06-25 [1] Github (r-spatial/stars@65c8726) #> styler 1.10.3 2024-04-07 [2] CRAN (R 4.4.0) #> units 0.8-5 2023-11-28 [2] CRAN (R 4.4.0) #> vctrs 0.6.5 2023-12-01 [2] CRAN (R 4.4.0) #> withr 3.0.0 2024-01-16 [2] CRAN (R 4.4.0) #> xfun 0.43 2024-03-25 [2] CRAN (R 4.4.0) #> yaml 2.3.8 2023-12-11 [2] CRAN (R 4.4.0) #> #> [1] C:/Users/500596972/AppData/Local/R/win-library/4.4 #> [2] C:/Program Files/R/R-4.4.0/library #> #> ------------------------------------------------------------------------------ ```