From b11bae1a65a5b9dcfdcb6c3b988038d242f5764d Mon Sep 17 00:00:00 2001 From: Sebastian Krantz Date: Fri, 7 Jul 2023 03:41:28 +0200 Subject: [PATCH] * Getting rid of dplyr and sp: replacing with collapse, sf and s2 + some optimization * Small fix. * Need to add sf and rename count column. * Using R's numeric() instead of alloc(). * Update NEWS * Fixing #63. The issue on M1 Mac is due to implicit conversion of a double (factor) to long int here: on M1 Macs this is rounded downwards, thus I needed to add a small number. Took me 6 freaking hours of simultaneous debugging on Mac and Windows to find it... * Adding return_sf argument: set to FALSE for memory efficient long datat frame (can be converted to S2 and other classes). * Appeasing R CMD check. * Comment and change of notation as requested by @r-barnes. * Using dgprintf as requested by @r-barnes. * Update CONTRIBUTING.md to explain C++ workflow * Also using dgprintf here. * Install instructions for fork. * Update README.md * Remove magrittr pipe. * Need to add dplyr for vignette it seems. --- CONTRIBUTING.md | 23 ++++++- DESCRIPTION | 27 +++----- NAMESPACE | 22 +++--- NEWS | 14 ++++ R/cwrapper.R | 140 +++++++++++++++++++------------------- R/dggridR.R | 154 ++++++++++++++++++++++++++---------------- R/zzz.R | 7 +- README.md | 17 ++++- man/dgcellstogrid.Rd | 4 +- man/dgearthgrid.Rd | 6 +- man/dgrectgrid.Rd | 7 +- man/dgshptogrid.Rd | 9 +-- src/DgDmdIDGG.cpp | 3 +- src/DgHexIDGG.cpp | 4 +- src/DgTriIDGG.cpp | 3 +- src/shputils.c | 4 +- vignettes/dggridR.Rmd | 17 +++-- 17 files changed, 273 insertions(+), 188 deletions(-) diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md index a4e41ac..5d86d65 100644 --- a/CONTRIBUTING.md +++ b/CONTRIBUTING.md @@ -11,4 +11,25 @@ These parts are: * `copy_to_src/func_gen.py` - Autogenerates full-factorial coordinate conversion code * `./update_from_upstream.h` - Copies the appropriate files from `submodules/DGGRID` into `src/` -Why do we have to do this copying? Because R's build model is pretty terrible: it assumes all the code lives in the same directory. If you want a reasonable, subdirectory-based organization of your code things become very painful. \ No newline at end of file +Why do we have to do this copying? Because R's build model is pretty terrible: it assumes all the code lives in the same directory. If you want a reasonable, subdirectory-based organization of your code things become very painful. + +Making changes to the C++ +============================================== + +Updating the C++ code in dggridR is a multi-step process. + +1. Pull code from https://github.com/sahrk/DGGRID into https://github.com/r-barnes/DGGRID +2. Update submodules/DGGRID from https://github.com/r-barnes/DGGRID +3. Run ./update_from_upstream.sh to update dggridR/src from submodules/DGGRID + +Making fixes to the C++ code requires running this process in reverse. The following workflow is suggested. + +1. Make any changes you need to submodules/DGGRID. Use submodules/DGGRID/CMakeLists.txt and the usual CMake workflow to ensure these compile. +2. Use ./update_from_upstream.sh to import your changes into dggridR and test that they work. +3. Iterate on the above until everything is working. +4. If you need changes to CRAN soonish, collect changes to submodules/DGGRID into a branch on r-barnes/DGGRID that submodules/DGGRID can point to. +5. Break changes in submodules/DGGRID apart into small branches on r-barnes/DGGRID which can be made into pull requests to https://github.com/sahrk/DGGRID. +6. Get Kevin Sahr's to approve and merge the PR's into https://github.com/sahrk/DGGRID. +7. Update https://github.com/r-barnes/DGGRID +8. Point submodules/DGGRID to master on https://github.com/r-barnes/DGGRID +9. Run ./update_from_upstream.sh to update dggridR/src from submodules/DGGRID \ No newline at end of file diff --git a/DESCRIPTION b/DESCRIPTION index 6f6dc56..d5242d0 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,19 +1,16 @@ Package: dggridR Type: Package Title: Discrete Global Grids -Version: 3.0.0 -Date: 2023-01-16 -Author: Richard Barnes [aut, cre], Kevin Sahr [aut, cph], Gerald Evenden [cph], Angus Johnson [cph], Frank Warmerdam [cph], Even Rouault [cph], Lian Song [ctb] +Version: 3.1.0 +Date: 2023-07-10 +Author: Richard Barnes [aut, cre], Kevin Sahr [aut, cph], Gerald Evenden [cph], Angus Johnson [cph], Frank Warmerdam [cph], Even Rouault [cph], Lian Song [ctb], Sebastian Krantz [ctb] Maintainer: Richard Barnes NeedsCompilation: yes -Depends: R (>= 3.4.0), dplyr (>= 0.4), rlang (>= 0.4), sf (>= 1.0), sp (>= 1.2) -Suggests: ggplot2, - knitr, - rmarkdown, - maps, - mapproj, - R.rsp, - testthat +Depends: R (>= 3.4.0) +Imports: Rcpp (>= 0.12.12), collapse (>= 1.8.0), sf (>= 1.0), s2 (>= 1.1) +LinkingTo: Rcpp +RcppModules: dgfuncs, gridgens, gridstats +Suggests: ggplot2, dplyr, knitr, rmarkdown, maps, mapproj, R.rsp, testthat License: AGPL (>= 3) Authors@R: c(person("Richard", "Barnes", role = c("aut", "cre"), email = "rijard.barnes@gmail.com"), person("Kevin", "Sahr", role = c("aut", "cph"), email = "sahrk@sou.edu"), @@ -21,16 +18,12 @@ Authors@R: c(person("Richard", "Barnes", role = c("aut", "cre"), email = "rij person("Angus", "Johnson", role = "cph"), person("Frank", "Warmerdam", role = "cph"), person("Even", "Rouault", role = "cph"), - person("Lian", "Song", role = "ctb") + person("Lian", "Song", role = "ctb"), + person("Sebastian", "Krantz", role = "ctb") ) Description: Spatial analyses involving binning require that every bin have the same area, but this is impossible using a rectangular grid laid over the Earth or over any projection of the Earth. Discrete global grids use hexagons, triangles, and diamonds to overcome this issue, overlaying the Earth with equally-sized bins. This package provides utilities for working with discrete global grids, along with utilities to aid in plotting such data. URL: https://github.com/r-barnes/dggridR/ BugReports: https://github.com/r-barnes/dggridR/ -Imports: - Rcpp (>= 0.12.12), - methods (>= 3.4.0) -LinkingTo: Rcpp -RcppModules: dgfuncs, gridgens, gridstats RoxygenNote: 7.2.3 Encoding: UTF-8 VignetteBuilder: knitr, R.rsp diff --git a/NAMESPACE b/NAMESPACE index 9427f42..cbc09c6 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -45,13 +45,17 @@ export(dgrectgrid) export(dgsetres) export(dgshptogrid) export(dgverify) -import(dplyr) -import(sf) -import(sp) -importFrom(methods,as) -importFrom(rlang,.data) -importFrom(utils,read.csv) -importFrom(utils,read.table) -importFrom(utils,tail) -importFrom(utils,write.table) +importFrom(Rcpp,loadModule) +importFrom(collapse,fgroup_by) +importFrom(collapse,fmutate) +importFrom(collapse,fsummarise) +importFrom(collapse,funique) +importFrom(collapse,qDF) +importFrom(s2,s2_convex_hull_agg) +importFrom(s2,s2_geog_point) +importFrom(sf,st_as_sf) +importFrom(sf,st_bbox) +importFrom(sf,st_read) +importFrom(sf,write_sf) +importFrom(tools,file_path_sans_ext) useDynLib(dggridR) diff --git a/NEWS b/NEWS index 5caf92a..feeb02b 100644 --- a/NEWS +++ b/NEWS @@ -1,3 +1,17 @@ +dggridR v3.1.0 +=============== +* Fixed a numeric precision issue on Apple silicon computers (#63) +* Removed dplyr, rlang and sp as dependencies to make the package more + lightweight and up to date (sp is depreciated) +* Added collapse and s2 as dependencies (incurring 0 additional dependencies) + yielding significant performance gains, especially in grid materialization +* Function dgshptogrid now also allows passing an sf data frame in memory, + in addition to a shapefile saved to disc +* Functions dgcellstogrid, dgearthgrid, dgrectgrid and dgshptogrid now have an + additional argument return_sf = TRUE. Setting return_sf = FALSE gives a more + memory efficient long data frame with the coordinates of the cell vertices. +* These changes were made by Sebastian Krantz + dggridR v3.0.0 (Release date: 2023-01-02) =============== * Fixes memory leaks in DGGRID diff --git a/R/cwrapper.R b/R/cwrapper.R index e9a4df5..f66aa4e 100644 --- a/R/cwrapper.R +++ b/R/cwrapper.R @@ -27,8 +27,8 @@ dgGEO_to_GEO <- function(dggs, in_lon_deg, in_lat_deg){ dgverify(dggs) N <- length(in_lon_deg) - out_lon_deg <- rep(0,N) - out_lat_deg <- rep(0,N) + out_lon_deg <- numeric(N) + out_lat_deg <- numeric(N) GEO_to_GEO(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_lon_deg, in_lat_deg, out_lon_deg, out_lat_deg) @@ -66,9 +66,9 @@ dgGEO_to_PROJTRI <- function(dggs, in_lon_deg, in_lat_deg){ dgverify(dggs) N <- length(in_lon_deg) - out_tnum <- rep(0,N) - out_tx <- rep(0,N) - out_ty <- rep(0,N) + out_tnum <- numeric(N) + out_tx <- numeric(N) + out_ty <- numeric(N) GEO_to_PROJTRI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_lon_deg, in_lat_deg, out_tnum, out_tx, out_ty) @@ -107,9 +107,9 @@ dgGEO_to_Q2DD <- function(dggs, in_lon_deg, in_lat_deg){ dgverify(dggs) N <- length(in_lon_deg) - out_quad <- rep(0,N) - out_qx <- rep(0,N) - out_qy <- rep(0,N) + out_quad <- numeric(N) + out_qx <- numeric(N) + out_qy <- numeric(N) GEO_to_Q2DD(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_lon_deg, in_lat_deg, out_quad, out_qx, out_qy) @@ -148,9 +148,9 @@ dgGEO_to_Q2DI <- function(dggs, in_lon_deg, in_lat_deg){ dgverify(dggs) N <- length(in_lon_deg) - out_quad <- rep(0,N) - out_i <- rep(0,N) - out_j <- rep(0,N) + out_quad <- numeric(N) + out_i <- numeric(N) + out_j <- numeric(N) GEO_to_Q2DI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_lon_deg, in_lat_deg, out_quad, out_i, out_j) @@ -189,7 +189,7 @@ dgGEO_to_SEQNUM <- function(dggs, in_lon_deg, in_lat_deg){ dgverify(dggs) N <- length(in_lon_deg) - out_seqnum <- rep(0,N) + out_seqnum <- numeric(N) GEO_to_SEQNUM(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_lon_deg, in_lat_deg, out_seqnum) @@ -226,8 +226,8 @@ dgGEO_to_PLANE <- function(dggs, in_lon_deg, in_lat_deg){ dgverify(dggs) N <- length(in_lon_deg) - out_px <- rep(0,N) - out_py <- rep(0,N) + out_px <- numeric(N) + out_py <- numeric(N) GEO_to_PLANE(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_lon_deg, in_lat_deg, out_px, out_py) @@ -266,8 +266,8 @@ dgPROJTRI_to_GEO <- function(dggs, in_tnum, in_tx, in_ty){ dgverify(dggs) N <- length(in_tnum) - out_lon_deg <- rep(0,N) - out_lat_deg <- rep(0,N) + out_lon_deg <- numeric(N) + out_lat_deg <- numeric(N) PROJTRI_to_GEO(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_tnum, in_tx, in_ty, out_lon_deg, out_lat_deg) @@ -306,9 +306,9 @@ dgPROJTRI_to_PROJTRI <- function(dggs, in_tnum, in_tx, in_ty){ dgverify(dggs) N <- length(in_tnum) - out_tnum <- rep(0,N) - out_tx <- rep(0,N) - out_ty <- rep(0,N) + out_tnum <- numeric(N) + out_tx <- numeric(N) + out_ty <- numeric(N) PROJTRI_to_PROJTRI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_tnum, in_tx, in_ty, out_tnum, out_tx, out_ty) @@ -348,9 +348,9 @@ dgPROJTRI_to_Q2DD <- function(dggs, in_tnum, in_tx, in_ty){ dgverify(dggs) N <- length(in_tnum) - out_quad <- rep(0,N) - out_qx <- rep(0,N) - out_qy <- rep(0,N) + out_quad <- numeric(N) + out_qx <- numeric(N) + out_qy <- numeric(N) PROJTRI_to_Q2DD(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_tnum, in_tx, in_ty, out_quad, out_qx, out_qy) @@ -390,9 +390,9 @@ dgPROJTRI_to_Q2DI <- function(dggs, in_tnum, in_tx, in_ty){ dgverify(dggs) N <- length(in_tnum) - out_quad <- rep(0,N) - out_i <- rep(0,N) - out_j <- rep(0,N) + out_quad <- numeric(N) + out_i <- numeric(N) + out_j <- numeric(N) PROJTRI_to_Q2DI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_tnum, in_tx, in_ty, out_quad, out_i, out_j) @@ -432,7 +432,7 @@ dgPROJTRI_to_SEQNUM <- function(dggs, in_tnum, in_tx, in_ty){ dgverify(dggs) N <- length(in_tnum) - out_seqnum <- rep(0,N) + out_seqnum <- numeric(N) PROJTRI_to_SEQNUM(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_tnum, in_tx, in_ty, out_seqnum) @@ -470,8 +470,8 @@ dgPROJTRI_to_PLANE <- function(dggs, in_tnum, in_tx, in_ty){ dgverify(dggs) N <- length(in_tnum) - out_px <- rep(0,N) - out_py <- rep(0,N) + out_px <- numeric(N) + out_py <- numeric(N) PROJTRI_to_PLANE(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_tnum, in_tx, in_ty, out_px, out_py) @@ -510,8 +510,8 @@ dgQ2DD_to_GEO <- function(dggs, in_quad, in_qx, in_qy){ dgverify(dggs) N <- length(in_quad) - out_lon_deg <- rep(0,N) - out_lat_deg <- rep(0,N) + out_lon_deg <- numeric(N) + out_lat_deg <- numeric(N) Q2DD_to_GEO(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_qx, in_qy, out_lon_deg, out_lat_deg) @@ -550,9 +550,9 @@ dgQ2DD_to_PROJTRI <- function(dggs, in_quad, in_qx, in_qy){ dgverify(dggs) N <- length(in_quad) - out_tnum <- rep(0,N) - out_tx <- rep(0,N) - out_ty <- rep(0,N) + out_tnum <- numeric(N) + out_tx <- numeric(N) + out_ty <- numeric(N) Q2DD_to_PROJTRI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_qx, in_qy, out_tnum, out_tx, out_ty) @@ -592,9 +592,9 @@ dgQ2DD_to_Q2DD <- function(dggs, in_quad, in_qx, in_qy){ dgverify(dggs) N <- length(in_quad) - out_quad <- rep(0,N) - out_qx <- rep(0,N) - out_qy <- rep(0,N) + out_quad <- numeric(N) + out_qx <- numeric(N) + out_qy <- numeric(N) Q2DD_to_Q2DD(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_qx, in_qy, out_quad, out_qx, out_qy) @@ -634,9 +634,9 @@ dgQ2DD_to_Q2DI <- function(dggs, in_quad, in_qx, in_qy){ dgverify(dggs) N <- length(in_quad) - out_quad <- rep(0,N) - out_i <- rep(0,N) - out_j <- rep(0,N) + out_quad <- numeric(N) + out_i <- numeric(N) + out_j <- numeric(N) Q2DD_to_Q2DI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_qx, in_qy, out_quad, out_i, out_j) @@ -676,7 +676,7 @@ dgQ2DD_to_SEQNUM <- function(dggs, in_quad, in_qx, in_qy){ dgverify(dggs) N <- length(in_quad) - out_seqnum <- rep(0,N) + out_seqnum <- numeric(N) Q2DD_to_SEQNUM(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_qx, in_qy, out_seqnum) @@ -714,8 +714,8 @@ dgQ2DD_to_PLANE <- function(dggs, in_quad, in_qx, in_qy){ dgverify(dggs) N <- length(in_quad) - out_px <- rep(0,N) - out_py <- rep(0,N) + out_px <- numeric(N) + out_py <- numeric(N) Q2DD_to_PLANE(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_qx, in_qy, out_px, out_py) @@ -754,8 +754,8 @@ dgQ2DI_to_GEO <- function(dggs, in_quad, in_i, in_j){ dgverify(dggs) N <- length(in_quad) - out_lon_deg <- rep(0,N) - out_lat_deg <- rep(0,N) + out_lon_deg <- numeric(N) + out_lat_deg <- numeric(N) Q2DI_to_GEO(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_i, in_j, out_lon_deg, out_lat_deg) @@ -794,9 +794,9 @@ dgQ2DI_to_PROJTRI <- function(dggs, in_quad, in_i, in_j){ dgverify(dggs) N <- length(in_quad) - out_tnum <- rep(0,N) - out_tx <- rep(0,N) - out_ty <- rep(0,N) + out_tnum <- numeric(N) + out_tx <- numeric(N) + out_ty <- numeric(N) Q2DI_to_PROJTRI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_i, in_j, out_tnum, out_tx, out_ty) @@ -836,9 +836,9 @@ dgQ2DI_to_Q2DD <- function(dggs, in_quad, in_i, in_j){ dgverify(dggs) N <- length(in_quad) - out_quad <- rep(0,N) - out_qx <- rep(0,N) - out_qy <- rep(0,N) + out_quad <- numeric(N) + out_qx <- numeric(N) + out_qy <- numeric(N) Q2DI_to_Q2DD(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_i, in_j, out_quad, out_qx, out_qy) @@ -878,9 +878,9 @@ dgQ2DI_to_Q2DI <- function(dggs, in_quad, in_i, in_j){ dgverify(dggs) N <- length(in_quad) - out_quad <- rep(0,N) - out_i <- rep(0,N) - out_j <- rep(0,N) + out_quad <- numeric(N) + out_i <- numeric(N) + out_j <- numeric(N) Q2DI_to_Q2DI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_i, in_j, out_quad, out_i, out_j) @@ -920,7 +920,7 @@ dgQ2DI_to_SEQNUM <- function(dggs, in_quad, in_i, in_j){ dgverify(dggs) N <- length(in_quad) - out_seqnum <- rep(0,N) + out_seqnum <- numeric(N) Q2DI_to_SEQNUM(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_i, in_j, out_seqnum) @@ -958,8 +958,8 @@ dgQ2DI_to_PLANE <- function(dggs, in_quad, in_i, in_j){ dgverify(dggs) N <- length(in_quad) - out_px <- rep(0,N) - out_py <- rep(0,N) + out_px <- numeric(N) + out_py <- numeric(N) Q2DI_to_PLANE(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_quad, in_i, in_j, out_px, out_py) @@ -996,8 +996,8 @@ dgSEQNUM_to_GEO <- function(dggs, in_seqnum){ dgverify(dggs) N <- length(in_seqnum) - out_lon_deg <- rep(0,N) - out_lat_deg <- rep(0,N) + out_lon_deg <- numeric(N) + out_lat_deg <- numeric(N) SEQNUM_to_GEO(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_seqnum, out_lon_deg, out_lat_deg) @@ -1034,9 +1034,9 @@ dgSEQNUM_to_PROJTRI <- function(dggs, in_seqnum){ dgverify(dggs) N <- length(in_seqnum) - out_tnum <- rep(0,N) - out_tx <- rep(0,N) - out_ty <- rep(0,N) + out_tnum <- numeric(N) + out_tx <- numeric(N) + out_ty <- numeric(N) SEQNUM_to_PROJTRI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_seqnum, out_tnum, out_tx, out_ty) @@ -1074,9 +1074,9 @@ dgSEQNUM_to_Q2DD <- function(dggs, in_seqnum){ dgverify(dggs) N <- length(in_seqnum) - out_quad <- rep(0,N) - out_qx <- rep(0,N) - out_qy <- rep(0,N) + out_quad <- numeric(N) + out_qx <- numeric(N) + out_qy <- numeric(N) SEQNUM_to_Q2DD(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_seqnum, out_quad, out_qx, out_qy) @@ -1114,9 +1114,9 @@ dgSEQNUM_to_Q2DI <- function(dggs, in_seqnum){ dgverify(dggs) N <- length(in_seqnum) - out_quad <- rep(0,N) - out_i <- rep(0,N) - out_j <- rep(0,N) + out_quad <- numeric(N) + out_i <- numeric(N) + out_j <- numeric(N) SEQNUM_to_Q2DI(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_seqnum, out_quad, out_i, out_j) @@ -1154,7 +1154,7 @@ dgSEQNUM_to_SEQNUM <- function(dggs, in_seqnum){ dgverify(dggs) N <- length(in_seqnum) - out_seqnum <- rep(0,N) + out_seqnum <- numeric(N) SEQNUM_to_SEQNUM(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_seqnum, out_seqnum) @@ -1190,8 +1190,8 @@ dgSEQNUM_to_PLANE <- function(dggs, in_seqnum){ dgverify(dggs) N <- length(in_seqnum) - out_px <- rep(0,N) - out_py <- rep(0,N) + out_px <- numeric(N) + out_py <- numeric(N) SEQNUM_to_PLANE(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], N, in_seqnum, out_px, out_py) diff --git a/R/dggridR.R b/R/dggridR.R index 969c0ba..ec03a83 100644 --- a/R/dggridR.R +++ b/R/dggridR.R @@ -1,10 +1,42 @@ -#' @import dplyr -#' @import sf -#' @import sp -#' @importFrom methods as -#' @importFrom rlang .data -#' @importFrom utils read.csv read.table tail write.table +#' @importFrom sf st_bbox st_as_sf write_sf st_read +#' @importFrom s2 s2_geog_point s2_convex_hull_agg +#' @importFrom collapse qDF fgroup_by fsummarise fmutate funique +#' @importFrom tools file_path_sans_ext #' @useDynLib dggridR +#' + +# Convert sf::st_bbox to sp::bbox +st_bbox_to_sp <- function(x) { + if(length(x) != 4L) stop("invalid bbox, needs to be a vector from st_bbox(), with exactly 4 elements") + matrix(x, ncol = 2, dimnames = list(c("x", "y"), c("min", "max"))) +} + +# Copied from sp to get rid of it as a dependency +makegrid <- function(bb, n = 10000, nsig = 2, cellsize, offset = rep(0.5, nrow(bb)), pretty = TRUE) { + # if (is(x, "Spatial")) + # bb = bbox(x) + # else bb = x + if (missing(cellsize)) { + pw = 1/nrow(bb) + cellsize = signif((prod(apply(bb, 1, diff))/n)^pw, nsig) + } + if (length(cellsize) == 1) cellsize = rep(cellsize, nrow(bb)) + min.coords = bb[, 1] + offset * cellsize + if (pretty) min.coords = signif(min.coords, max(ceiling(log10(abs(bb[1, ])/cellsize)))) + sel = min.coords - offset * cellsize > bb[, 1] + if (any(sel)) min.coords[sel] = min.coords[sel] - cellsize[sel] + expand.grid.arglist = list() + for (i in 1:nrow(bb)) { + name = paste("x", i, sep = "") + from = min.coords[i] + by = cellsize[i] + length.out = round(1 + (bb[i, 2] - from)/by) + expand.grid.arglist[[name]] = seq(from, by = by, length.out = length.out) + } + xy = do.call(expand.grid, expand.grid.arglist) + attr(xy, "cellsize") = cellsize + return(xy) +} #' @name dg_env @@ -504,14 +536,24 @@ dg_closest_res_to_cls <- function(dggs,cls,round='nearest',show_info=TRUE,metric #' #' @return Returns an sf object. #' -dg_process_polydata <- function(polydata){ - polydata <- as.data.frame(polydata) - - polydata %>% - sf::st_as_sf(coords=c("x", "y"), crs=4326) %>% - group_by(.data$seqnum) %>% - summarise(geometry = sf::st_combine(geometry)) %>% - sf::st_cast("POLYGON") +dg_process_polydata <- function(polydata) { + + x <- y <- seqnum <- geometry <- NULL # For R CMD Check: no visible binding for global variables + + # Using s2 directly: faster !! + qDF(polydata) |> + fmutate(geometry = s2_geog_point(x, y)) |> + fgroup_by(seqnum, sort = TRUE) |> + fsummarise(geometry = s2_convex_hull_agg(geometry)) |> + st_as_sf(crs = 4326) + + # sf solution: also much faster than with dplyr, but 3x slower than s2 + # qDF(polydata) |> + # st_as_sf(coords = c("x", "y"), crs = 4326) |> + # fgroup_by(seqnum, sort = TRUE) |> + # fsummarise(geometry = st_combine(geometry)) |> + # fmutate(geometry = `oldClass<-`(geometry, c("sfc_MULTIPOINT", "sfc"))) |> + # st_cast("POLYGON") } @@ -539,10 +581,7 @@ dg_process_polydata <- function(polydata){ #' generate the grid. Small values yield long generation times #' while large values may omit cells. #' -#' @param savegrid If savegrid is set to a file path, then a shapefile -#' containing the grid is written to that path and the filename -#' is returned. No other manipulations are done. -#' Default: NA (do not save grid, return it) +#' @param \dots Further arguments passed to \code{\link{dgcellstogrid}}. #' #' @return Returns an sf object. #' If \code{!is.na(savegrid)}, returns a filename. @@ -556,20 +595,18 @@ dg_process_polydata <- function(polydata){ #' minlat=24.7433195, minlon=-124.7844079, #' maxlat=49.3457868, maxlon=-66.9513812) #' @export -dgrectgrid <- function(dggs,minlat=-1,minlon=-1,maxlat=-1,maxlon=-1,cellsize=0.1,savegrid=NA){ #TODO: Densify? +dgrectgrid <- function(dggs,minlat=-1,minlon=-1,maxlat=-1,maxlon=-1,cellsize=0.1, ...){ #TODO: Densify? dgverify(dggs) - coords <- matrix(c(minlon, minlat, maxlon, minlat, maxlon, maxlat, maxlon, minlat, minlon, minlat), ncol = 2, byrow = TRUE) - regbox <- Polygon(coords) - regbox <- SpatialPolygons(list(Polygons(list(regbox), ID = "a")), proj4string=CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")) + bbox_coords <- matrix(c(minlon, minlat, maxlon, maxlat), ncol = 2, dimnames = list(c("x", "y"), c("min", "max"))) #Generate a dense grid of points - samp_points <- sp::makegrid(regbox, cellsize = cellsize) + samp_points <- makegrid(bbox_coords, cellsize = cellsize) #Convert the points to SEQNUM ids for dggs - samp_points <- dgGEO_to_SEQNUM(dggs,samp_points$x1, samp_points$x2)$seqnum + samp_points <- dgGEO_to_SEQNUM(dggs, samp_points$x1, samp_points$x2)$seqnum - dgcellstogrid(dggs, samp_points, savegrid=savegrid) + dgcellstogrid(dggs, samp_points, ...) } @@ -581,13 +618,11 @@ dgrectgrid <- function(dggs,minlat=-1,minlon=-1,maxlat=-1,maxlon=-1,cellsize=0.1 #' #' @description Note: If you have a high-resolution grid this may take a #' very long time to execute. +#' +#' +#' @param dggs A dggs object from dgconstruct(). +#' @inheritParams dgcellstogrid #' -#' @param dggs A dggs object from dgconstruct() -#' -#' @param savegrid If savegrid is set to a file path, then a shapefile -#' containing the grid is written to that path and the filename -#' is returned. No other manipulations are done. -#' Default: NA (do not save grid, return it) #' #' @return Returns an sf object. #' If \code{!is.na(savegrid)}, returns a filename. @@ -601,11 +636,12 @@ dgrectgrid <- function(dggs,minlat=-1,minlon=-1,maxlat=-1,maxlon=-1,cellsize=0.1 #' gridfilename <- dgearthgrid(dggs,savegrid=tempfile(fileext=".shp")) #Save directly to a file #' } #' @export -dgearthgrid <- function(dggs,savegrid=NA){ #TODO: Densify? +dgearthgrid <- function(dggs, savegrid = NA, return_sf = TRUE) { #TODO: Densify? dgverify(dggs) grid <- GlobalGrid(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]]) - if(is.na(savegrid)){ + if(is.na(savegrid)) { + if(!return_sf) return(qDF(grid)) dg_process_polydata(grid) } else { grid <- dg_process_polydata(grid) @@ -631,6 +667,8 @@ dgearthgrid <- function(dggs,savegrid=NA){ #TODO: Densify? #' containing the grid is written to that path and the filename #' is returned. No other manipulations are done. #' Default: NA (do not save grid, return it) +#' +#' @param return_sf logical. If \code{FALSE}, a long-format data frame giving the coordinates of the vertices of each cell is returned. This is is considerably faster and more memory efficient than creating an sf data frame. #' #' @return Returns an sf object. #' If \code{!is.na(savegrid)}, returns a filename. @@ -646,19 +684,20 @@ dgearthgrid <- function(dggs,savegrid=NA){ #TODO: Densify? #' #Get grid cells for the earthquakes identified #' grid <- dgcellstogrid(dggs, dgquakes$cell) #' @export -dgcellstogrid <- function(dggs, cells, savegrid=NA){ #TODO: Densify? +dgcellstogrid <- function(dggs, cells, savegrid=NA, return_sf = TRUE){ #TODO: Densify? dgverify(dggs) #dggrid also eliminates duplicate cells, but doing so here saves disk space #and likely wall time, given the costs of IO, not that it matters unless the #data set is huge - cells <- unique(cells) + cells <- funique(cells) - if(max(cells)>dgmaxcell(dggs)) + if(max(cells) > dgmaxcell(dggs)) stop("'cells' contained cell ids which were larger than the maximum id!") grid <- SeqNumGrid(dggs[["pole_lon_deg"]], dggs[["pole_lat_deg"]], dggs[["azimuth_deg"]], dggs[["aperture"]], dggs[["res"]], dggs[["topology"]], dggs[["projection"]], cells) if(is.na(savegrid)){ + if(!return_sf) return(qDF(grid)) dg_process_polydata(grid) } else { grid <- dg_process_polydata(grid) @@ -679,7 +718,7 @@ dgcellstogrid <- function(dggs, cells, savegrid=NA){ #TODO: Densify? #' #' @return The filename the grid was saved to dgsavegrid <- function(grid,shpfname) { - sf::write_sf(grid, shpfname, driver='ESRI Shapefile', layer='dggrid') + write_sf(grid, shpfname, driver='ESRI Shapefile', layer='dggrid') shpfname } @@ -710,16 +749,13 @@ dgsavegrid <- function(grid,shpfname) { #' #' @param dggs A dggs object from dgconstruct() #' -#' @param shpfname File name of the shapefile. Filename should end with '.shp' +#' @param shpfname Either a sf data frame or the file name of the shapefile. Filename should end with '.shp'. #' #' @param cellsize Distance, in degrees, between the sample points used to #' generate the grid. Small values yield long generation times #' while large values may omit cells. #' -#' @param savegrid If savegrid is set to a file path, then a shapefile -#' containing the grid is written to that path and the filename -#' is returned. No other manipulations are done. -#' Default: NA (do not save grid, return it) +#' @param \dots Further arguments passed to \code{\link{dgcellstogrid}}. #' #' @return Returns an sf object. #' If \code{!is.na(savegrid)}, returns a filename. @@ -730,25 +766,29 @@ dgsavegrid <- function(grid,shpfname) { #' dggs <- dgconstruct(spacing=25, metric=FALSE, resround='nearest') #' south_africa_grid <- dgshptogrid(dggs,dg_shpfname_south_africa()) #' @export -dgshptogrid <- function(dggs,shpfname,cellsize=0.1,savegrid=NA){ #TODO: Densify? +dgshptogrid <- function(dggs, shpfname, cellsize = 0.1, ...) { #TODO: Densify? dgverify(dggs) - - shpfname <- trimws(shpfname) - - if(!grepl('\\.shp$',shpfname)) - stop("Shapefile name does to end with '.shp'!") - if(!file.exists(shpfname)) - stop('Shapefile does not exist!') - - dsn <- dirname(shpfname) - layer <- tools::file_path_sans_ext(basename(shpfname)) - poly <- sf::st_read(dsn, layer=layer) + + if(!(is.data.frame(shpfname) && inherits(shpfname, "sf"))) { + + shpfname <- trimws(shpfname) + + if(!grepl('\\.shp$',shpfname)) + stop("Shapefile name does to end with '.shp'!") + if(!file.exists(shpfname)) + stop('Shapefile does not exist!') + + dsn <- dirname(shpfname) + layer <- file_path_sans_ext(basename(shpfname)) + bb <- st_bbox_to_sp(st_bbox(st_read(dsn, layer=layer))) + + } else bb <- st_bbox_to_sp(st_bbox(shpfname)) #Generate a dense grid of points - samp_points <- sp::makegrid(sf::as_Spatial(poly), cellsize = cellsize) + samp_points <- makegrid(bb, cellsize = cellsize) #Convert the points to SEQNUM ids for dggs - samp_points <- dgGEO_to_SEQNUM(dggs,samp_points$x1, samp_points$x2)$seqnum + samp_points <- dgGEO_to_SEQNUM(dggs, samp_points$x1, samp_points$x2)$seqnum - dgcellstogrid(dggs, samp_points, savegrid=savegrid) + dgcellstogrid(dggs, samp_points, ...) } diff --git a/R/zzz.R b/R/zzz.R index 1809752..d6edd3f 100644 --- a/R/zzz.R +++ b/R/zzz.R @@ -1,6 +1,7 @@ -Rcpp::loadModule("dgfuncs", TRUE) -Rcpp::loadModule("gridgens", TRUE) -Rcpp::loadModule("gridstats", TRUE) +#' @importFrom Rcpp loadModule +loadModule("dgfuncs", TRUE) +loadModule("gridgens", TRUE) +loadModule("gridstats", TRUE) #TODO: Used to display start-up messages #.onAttach <- function(...) { diff --git a/README.md b/README.md index 38236ee..d21517c 100644 --- a/README.md +++ b/README.md @@ -2,6 +2,19 @@ dggridR: Discrete Global Grids for R ==================================== +*** +*Note*: This fork implements vital fixes to *dggridR*, and performance improvements. The upstream project is currently unmaintained. This fork can be installed using either of: + +```r +remotes::install_github("SebKrantz/dggridR") + +install.packages("dggridR", + repos = c(ropensci = 'https://fastverse.r-universe.dev', + CRAN = 'https://cloud.r-project.org')) +``` + + +*** _Spatial Analysis Done Right_ @@ -55,7 +68,7 @@ Your analysis could be as easy as this: ```R library(dggridR) -library(dplyr) +library(collapse) #Construct a global grid with cells approximately 1000 miles across dggs <- dgconstruct(spacing=1000, metric=FALSE, resround='down') @@ -67,7 +80,7 @@ data(dgquakes) dgquakes$cell <- dgGEO_to_SEQNUM(dggs, dgquakes$lon, dgquakes$lat)$seqnum #Get the number of earthquakes in each equally-sized cell -quakecounts <- dgquakes %>% group_by(cell) %>% summarise(count=n()) +quakecounts <- dgquakes %>% fcount(cell) ``` Show me more examples! diff --git a/man/dgcellstogrid.Rd b/man/dgcellstogrid.Rd index 0f20f98..5380074 100644 --- a/man/dgcellstogrid.Rd +++ b/man/dgcellstogrid.Rd @@ -4,7 +4,7 @@ \alias{dgcellstogrid} \title{Return boundary coordinates for specified cells} \usage{ -dgcellstogrid(dggs, cells, savegrid = NA) +dgcellstogrid(dggs, cells, savegrid = NA, return_sf = TRUE) } \arguments{ \item{dggs}{A dggs object from dgconstruct()} @@ -15,6 +15,8 @@ dgcellstogrid(dggs, cells, savegrid = NA) containing the grid is written to that path and the filename is returned. No other manipulations are done. Default: NA (do not save grid, return it)} + +\item{return_sf}{logical. If \code{FALSE}, a long-format data frame giving the coordinates of the vertices of each cell is returned. This is is considerably faster and more memory efficient than creating an sf data frame.} } \value{ Returns an sf object. diff --git a/man/dgearthgrid.Rd b/man/dgearthgrid.Rd index f7b7e44..31e4e7a 100644 --- a/man/dgearthgrid.Rd +++ b/man/dgearthgrid.Rd @@ -5,15 +5,17 @@ \title{Return the coordinates constituting the boundary of cells for the entire Earth} \usage{ -dgearthgrid(dggs, savegrid = NA) +dgearthgrid(dggs, savegrid = NA, return_sf = TRUE) } \arguments{ -\item{dggs}{A dggs object from dgconstruct()} +\item{dggs}{A dggs object from dgconstruct().} \item{savegrid}{If savegrid is set to a file path, then a shapefile containing the grid is written to that path and the filename is returned. No other manipulations are done. Default: NA (do not save grid, return it)} + +\item{return_sf}{logical. If \code{FALSE}, a long-format data frame giving the coordinates of the vertices of each cell is returned. This is is considerably faster and more memory efficient than creating an sf data frame.} } \value{ Returns an sf object. diff --git a/man/dgrectgrid.Rd b/man/dgrectgrid.Rd index 8010f1e..e653c58 100644 --- a/man/dgrectgrid.Rd +++ b/man/dgrectgrid.Rd @@ -12,7 +12,7 @@ dgrectgrid( maxlat = -1, maxlon = -1, cellsize = 0.1, - savegrid = NA + ... ) } \arguments{ @@ -30,10 +30,7 @@ dgrectgrid( generate the grid. Small values yield long generation times while large values may omit cells.} -\item{savegrid}{If savegrid is set to a file path, then a shapefile -containing the grid is written to that path and the filename -is returned. No other manipulations are done. -Default: NA (do not save grid, return it)} +\item{\dots}{Further arguments passed to \code{\link{dgcellstogrid}}.} } \value{ Returns an sf object. diff --git a/man/dgshptogrid.Rd b/man/dgshptogrid.Rd index 72d2eca..c34e8fe 100644 --- a/man/dgshptogrid.Rd +++ b/man/dgshptogrid.Rd @@ -5,21 +5,18 @@ \title{Return boundary coordinates for cells intersecting a shapefile} \usage{ -dgshptogrid(dggs, shpfname, cellsize = 0.1, savegrid = NA) +dgshptogrid(dggs, shpfname, cellsize = 0.1, ...) } \arguments{ \item{dggs}{A dggs object from dgconstruct()} -\item{shpfname}{File name of the shapefile. Filename should end with '.shp'} +\item{shpfname}{Either a sf data frame or the file name of the shapefile. Filename should end with '.shp'.} \item{cellsize}{Distance, in degrees, between the sample points used to generate the grid. Small values yield long generation times while large values may omit cells.} -\item{savegrid}{If savegrid is set to a file path, then a shapefile -containing the grid is written to that path and the filename -is returned. No other manipulations are done. -Default: NA (do not save grid, return it)} +\item{\dots}{Further arguments passed to \code{\link{dgcellstogrid}}.} } \value{ Returns an sf object. diff --git a/src/DgDmdIDGG.cpp b/src/DgDmdIDGG.cpp index e788320..2273eaa 100644 --- a/src/DgDmdIDGG.cpp +++ b/src/DgDmdIDGG.cpp @@ -123,7 +123,8 @@ DgDmdIDGG::initialize (void) double factor = parentScaleFac * 2.0L; // aperture 4 scaleFac_ = factor; - maxD_ = factor - 1.0L; + // Adding small number (1e-6) to prevent rounding down in conversion to integer (fixes issue #63 experienced on Apple ARM computers) + maxD_ = factor+1e-6 - 1.0L; //cout << res() << " " << aperture(); //cout << " f: " << factor << " maxD: " << maxD_ << endl; diff --git a/src/DgHexIDGG.cpp b/src/DgHexIDGG.cpp index 463a507..a228dec 100644 --- a/src/DgHexIDGG.cpp +++ b/src/DgHexIDGG.cpp @@ -166,8 +166,8 @@ DgHexIDGG::initialize (void) if (isClassIII()) factor *= M_SQRT7; - - maxD_ = factor - 1.0; + // Adding small number (1e-6) to prevent rounding down in conversion to integer (fixes issue #63 experienced on Apple ARM computers) + maxD_ = factor+1e-6 - 1.0; //cout << res() << " " << aperture(); //cout << " f: " << factor << " maxD: " << maxD_ << endl; diff --git a/src/DgTriIDGG.cpp b/src/DgTriIDGG.cpp index 7a855ef..f5869ca 100644 --- a/src/DgTriIDGG.cpp +++ b/src/DgTriIDGG.cpp @@ -117,7 +117,8 @@ DgTriIDGG::initialize (void) double factor = parentScaleFac * 2.0L; // aperture 4 scaleFac_ = factor; - maxD_ = factor - 1.0L; + // Adding small number (1e-6) to prevent rounding down in conversion to integer (fixes issue #63 experienced on Apple ARM computers) + maxD_ = factor+1e-6 - 1.0L; //cout << res() << " " << aperture(); //cout << " f: " << factor << " maxD: " << maxD_ << endl; diff --git a/src/shputils.c b/src/shputils.c index fc80997..088a970 100644 --- a/src/shputils.c +++ b/src/shputils.c @@ -441,10 +441,10 @@ void showitems(void) dsum = dsum + dtmp; } mean=dsum/maxrec; - dgprintf(stmp,"%%.%df to %%.%df \t(%%.%df)",iDecimals,iDecimals,iDecimals); + dgprintf(stmp, "%%.%df to %%.%df \t(%%.%df)",iDecimals,iDecimals,iDecimals); if (dlow < dhigh) dgprintf(stmp,dlow,dhigh,mean); else if (dlow == dhigh) { - dgprintf(stmp,"= %%.%df",iDecimals); + dgprintf(stmp, "= %%.%df",iDecimals); dgprintf(stmp,dlow); } else dgprintf("No Values"); diff --git a/vignettes/dggridR.Rmd b/vignettes/dggridR.Rmd index e022a93..3e3fcb0 100644 --- a/vignettes/dggridR.Rmd +++ b/vignettes/dggridR.Rmd @@ -15,6 +15,7 @@ vignette: > ```{r, fig.width=5, fig.height=5, results='hide', warning=FALSE, error=FALSE, message=FALSE, echo=FALSE, fig.align='center'} # Generate cover picture library(dggridR) +library(sf) library(ggplot2) # Generate grids of various sizes @@ -26,11 +27,11 @@ countries <- map_data("world") # Crop generate dgrids to areas of interest bounds = st_bbox(c(xmin = -90, xmax = 75, ymin = -90, ymax = 90), crs = st_crs(4326)) -hgrids[[1]] = hgrids[[1]] %>% st_make_valid() %>% st_filter(st_as_sfc(bounds), .predicate=st_within) +hgrids[[1]] = hgrids[[1]] |> st_make_valid() |> st_filter(st_as_sfc(bounds), .predicate=st_within) bounds = st_bbox(c(xmin = 20, xmax = 145, ymin = -90, ymax = 90), crs = st_crs(4326)) -hgrids[[2]] = hgrids[[2]] %>% st_make_valid() %>% st_filter(st_as_sfc(bounds), .predicate=st_within) +hgrids[[2]] = hgrids[[2]] |> st_make_valid() |> st_filter(st_as_sfc(bounds), .predicate=st_within) bounds = st_bbox(c(xmin = 90, xmax = 215, ymin = -90, ymax = 90), crs = st_crs(4326)) -hgrids[[3]] = hgrids[[3]] %>% st_make_valid() %>% st_filter(st_as_sfc(bounds), .predicate=st_within) +hgrids[[3]] = hgrids[[3]] |> st_make_valid() |> st_filter(st_as_sfc(bounds), .predicate=st_within) ggplot() + geom_polygon(data=countries, aes(x=long, y=lat, group=group), fill=NA, color="black") + @@ -192,7 +193,7 @@ example demonstrates how to get the center coordinates of the cells. ```{r, results='hide', warning=FALSE, error=FALSE, message=FALSE} #Include libraries library(dggridR) -library(dplyr) +library(collapse) #Construct a global grid with cells approximately 1000 miles across dggs <- dgconstruct(spacing=1000, metric=FALSE, resround='down') @@ -207,7 +208,7 @@ dgquakes$cell <- dgGEO_to_SEQNUM(dggs,dgquakes$lon,dgquakes$lat)$seqnum cellcenters <- dgSEQNUM_to_GEO(dggs,dgquakes$cell) #Get the number of earthquakes in each cell -quakecounts <- dgquakes %>% group_by(cell) %>% summarise(count=n()) +quakecounts <- dgquakes |> fcount(cell, name = "count") #Get the grid cell boundaries for cells which had quakes grid <- dgcellstogrid(dggs,quakecounts$cell) @@ -218,8 +219,8 @@ grid <- merge(grid,quakecounts,by.x="seqnum",by.y="cell") #Make adjustments so the output is more visually interesting grid$count <- log(grid$count) -cutoff <- quantile(grid$count,0.9) -grid <- grid %>% mutate(count=ifelse(count>cutoff,cutoff,count)) +cutoff <- fquantile(grid$count, 0.9) +grid <- grid |> fmutate(count = ifelse(count>cutoff,cutoff,count)) #Get polygons for each country of the world countries <- map_data("world") @@ -297,7 +298,6 @@ their associated grid cells: ```{r, results='hide', warning=FALSE, error=FALSE, message=FALSE} #Include libraries library(dggridR) -library(dplyr) N <- 100 #How many cells to sample @@ -367,7 +367,6 @@ and generate a grid accordingly. ```{r, results='hide', warning=FALSE, error=FALSE, message=FALSE} #Include libraries library(dggridR) -library(dplyr) N <- 100 #How many cells to sample