-
Notifications
You must be signed in to change notification settings - Fork 3
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
[WIP] setting up Inst/Include header #18
Conversation
Oh sweet thanks so much. Hoping to level up and contribute more soonish! |
I've set up the header, and made Rcpp wrappers for the I've written a proof-of-concept x <- c(0, 0, 0.75, 1, 0.5, 0.8, 0.69)
y <- c(0, 1, 1, 0.8, 0.7, 0.6, 0)
(ind <- earcut(cbind(x, y)))
# [1] 2 1 7 7 6 5 5 4 3 2 7 5 5 3 2
plot_ears(cbind(x, y), ind)
mat <- cbind(x,y)
sfg <- sfheaders::sfg_polygon( mat, close = F )
(ind <- decido:::earcut_sfc( sfg )) + 1
# [1] 2 1 7 7 6 5 5 4 3 2 7 5 5 3 2
plot_ears(cbind(x, y), ind + 1) |
Here's a proof-of-concept on how it could be called by library(Rcpp)
cppFunction(
depends = c("decido","sfheaders","geometries") ## geometries because I'm using dev version of sfheaders
, includes = c(
'#include "decido/decido.hpp"'
, '#include "sfheaders/df/sfg.hpp"'
, '#include "sfheaders/sfheaders.hpp"'
)
, code = '
SEXP triangles( SEXP sfg ) {
Polygons polyrings = Rcpp::as< Polygons >( sfg );
std::vector< uint32_t > indices = mapbox::earcut< uint32_t >( polyrings );
Rcpp::IntegerVector idx = Rcpp::wrap( indices );
Rcpp::DataFrame df = sfheaders::df::sfg_to_df( sfg );
Rcpp::NumericVector x = df["x"];
Rcpp::NumericVector y = df["y"];
Rcpp::NumericVector xx = x[ idx ];
Rcpp::NumericVector yy = y[ idx ];
int n_triangles = idx.length() / 3;
Rcpp::IntegerVector triangle_idx = Rcpp::seq(1, n_triangles );
Rcpp::IntegerVector triangle_ids = Rcpp::rep_each( triangle_idx, 3 );
Rcpp::DataFrame df_tri = Rcpp::DataFrame::create(
Rcpp::_["triangle_id"] = triangle_ids,
Rcpp::_["x"] = xx,
Rcpp::_["y"] = yy
);
Rcpp::StringVector geometry_cols {"x","y"};
Rcpp::String polygon_id = "triangle_id";
std::string xyzm = "XY";
bool close = false;
Rcpp::List sfg_tri = sfheaders::sfc::sfc_polygon(
df_tri, geometry_cols, polygon_id, polygon_id, xyzm, close
);
return sfg_tri;
}
'
)
library(sf)
library(mapdeck)
set_token( secret::get_secret("mapbox"))
nc <- sf::st_read(system.file("./shape/nc.shp", package = "sf"))
sf_poly <- sfheaders::sf_cast( nc, "POLYGON" )
sfg <- sf_poly[1, ]$geometry[[1]]
sf_tri <- triangles( sfg )
sf <- sf::st_as_sf( sf_tri )
sf$id <- 1:nrow( sf )
mapdeck() %>%
add_polygon(
data = sf
, fill_colour = "id"
, stroke_width = 20
, tooltip = "id"
, legend = T
)
|
Note to me for future reference. I tried making a sf <- geojsonsf::geojson_sf("https://symbolixau.github.io/data/geojson/SA2_2016_VIC.json")
sf_poly <- sfheaders::sf_cast( sf, "POLYGON" )
sfc <- sf_poly$geometry
lst <- decido:::earcut_sfc( sfc )
lst2 <- decido:::earcut_sfc2( sfc )
all.equal(lst, lst2)
# [1] TRUE
microbenchmark::microbenchmark(
wrap = { lst <- decido:::earcut_sfc( sfc ) },
rcpp = { lst2 <- decido:::earcut_sfc2( sfc ) }
)
# Unit: milliseconds
# expr min lq mean median uq max neval
# wrap 203.1208 210.5301 225.9748 213.9351 218.7850 311.4982 100
# rcpp 285.7782 305.3584 337.2406 311.4472 385.1729 473.1197 100 (on branch inherited) |
@mdsumner I think with this last commit I've finished what I needed to do so is ready for your review. Things to note
Examplessf <- geojsonsf::geojson_sf("https://symbolixau.github.io/data/geojson/SA2_2016_VIC.json")
sf_poly <- sfheaders::sf_cast( sf, "POLYGON" )
library(Rcpp)
## Single Polygon
cppFunction(
depends = "decido"
, includes = '#include "decido/decido.hpp"'
, code = '
Rcpp::IntegerVector earcut( Rcpp::List polygon ) {
return decido::api::earcut( polygon );
}
'
)
## Single Polygon
which( sapply( sf_poly$geometry, length ) > 1 )
r <- 17
idx <- earcut( sf_poly[r, ]$geometry[[1]] )
df <- sfheaders::sf_to_df( sf[r, ])
ids <- rep(1:(length(idx)/3), each = 3)
df <- df[ idx, ]
df$id <- ids
sf_triangles <- sfheaders::sf_polygon(
obj = df
, polygon_id = "id"
, x = "x"
, y = "y"
)
mapdeck() %>%
add_polygon(
data = sf_triangles
, fill_colour = "id"
, fill_opacity = 0.6
, tooltip = "id"
) ### Multiple Polygons
cppFunction(
depends = "decido"
, includes = '#include "decido/decido.hpp"'
, code = '
Rcpp::List earcuts( Rcpp::List polygons ) {
R_xlen_t i;
R_xlen_t n = polygons.length();
Rcpp::List res( n );
for( i = 0; i < n; ++i ) {
Rcpp::List polygon = polygons[ i ];
res[i] = decido::api::earcut( polygon );
}
return res;
}
'
)
lst <- earcuts( sf_poly$geometry )
str( lst )
# List of 616
# $ : int [1:270] 92 91 90 90 89 88 88 87 86 86 ...
# $ : int [1:669] 225 224 223 223 222 221 220 219 218 218 ...
# $ : int [1:999] 333 332 331 331 330 329 329 328 327 327 ...
# $ : int [1:858] 1 288 287 287 286 285 285 284 283 281 ...
# $ : int [1:1626] 544 543 542 542 541 540 539 538 537 537 ...
# $ : int [1:258] 1 88 87 87 86 85 85 84 83 83 ...
|
Well, great day's work (!) - I can't see any problems. Thanks! I definitely think it's worth exposing that new function. Do you want to suggest a NEWS item? I did have one drafted earlier today but you might as well write it. It's Saturday night and I'm unlikely to delve much but really excited with what you've done here. I'll merge now if you want, it's passing check and hasn't changed anything that I was doing before afaics. |
yeah me either now I'm two glasses of wine down. But I'll write the NEWs entry tomorrow. |
yeah let's step away haha |
I've updated the NEWS and added some tests. I think the task of "exposing to R" is better suited to the discussion in #9 , as I think it needs a philosophical discussion on what the inputs and outputs are going to be, and which other (if any) packages are required as dependencies. Therefore, I'm inclined to say this PR is complete. |
Brilliant, thanks so much! No excuses for me now I have to get into this headers stuff :) |
I'm making a work-in-progress PR so you can see what I'm doing as and when I commit changes.
I'll let you know when I think it's ready
TODO
LinkingTo
from Rcpp code / packagesexposeearcut( Rcpp::List polygon )