Skip to content
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

sf prototype #9

Closed
mdsumner opened this issue Mar 15, 2018 · 11 comments
Closed

sf prototype #9

mdsumner opened this issue Mar 15, 2018 · 11 comments

Comments

@mdsumner
Copy link
Member

mdsumner commented Mar 15, 2018

Bare bones earcut for sf, for use in another package

#' Triangulate simple features
#' 
#' Ear clipping triangulation applied to simple features POLYGON. 
#'
#' @param x simple features object
#' @param ... ignored
#' @examples 
#' library(sf)
#' library(decido)
#' example(st_read)
#' earcut(nc)
#' earcut(st_cast(nc[1:10, ], "POLYGON"))
earcut.sf <- function(xy, ...) {
   geom <- xy[[attr(xy, "sf_column")]]
  if (!(inherits(geom, "sfc_POLYGON") || inherits(geom, "sfc_MULTIPOLYGON"))) stop("non POLYGON type")
  sf::st_set_geometry(xy, earcut(geom))
}

earcut.sfc <- function(xy, ...) {
  out <- lapply(xy, earcut)
##  sf::st_geometrycollection(lapply(out, function(ix) sf::st_polygon(list(ix))))
  sf::st_sfc(lapply(out, function(feature) sf::st_geometrycollection(lapply(feature, function(triangle) sf::st_polygon(list(triangle))))), 
             crs = sf::st_crs(xy))
}

earcut.POLYGON <- function(xy, ...) {
  ni <- unlist(lapply(xy, nrow))
  h <- head(cumsum(ni), -1) + 1
  if (length(h) < 1) h <- 0
  xy <- do.call(rbind, xy)[,1:2]  ## only 2D supported
  idx <- earcut(xy, h)
  lapply(split(idx, rep(seq_len(length(idx)/3), each = 3)), function(i) xy[c(i, i[1]), ])
 
}

earcut.MULTIPOLYGON <- function(xy, ...) {
  unlist(lapply(xy, earcut.POLYGON), recursive = FALSE)
  
}
@mdsumner mdsumner added the help wanted Extra attention is needed label Mar 15, 2018
@mdsumner
Copy link
Member Author

I tested this pretty heavily in the 0.2.0 checks, it definitely works!

dcooley added a commit to dcooley/decido that referenced this issue May 3, 2020
@dcooley
Copy link
Contributor

dcooley commented May 3, 2020

I have another work-in-progress in-bound for this issue

Things to consider

  • should it link to sfheaders, or do you want a stand-alone implementation so decido doesn't import or rely on any other packages?
  • I haven't benchmarked my code for speed. I'm sure there are optimisations I can make.
  • What object should it reutrn? I've currently got it returning MULTIPOLYGONS.
  • Will need to work on sf to so all the data attributes can be kept / replicated for each triangle
// [[Rcpp::export]]
SEXP earcut_sfc( Rcpp::List& sfc ) {
  // expecting an sfc objet

  Rcpp::List attributes = sfheaders::sfc::get_sfc_attributes( sfc );

  // TODO
  // - cast to POLYGON
  // - for each sfg
  // - tringulate
  // - convert to sfg_POLYGON
  R_xlen_t i;
  R_xlen_t n = sfc.length();
  Rcpp::List res( n );
  for( i = 0; i < n; ++i ) {
    SEXP poly = sfc[ i ];
    Rcpp::DataFrame df = sfheaders::df::sfg_to_df( poly );
    Rcpp::IntegerVector indices = decido::api::earcut( poly );
    R_xlen_t n_triangles = indices.length() / 3;

    Rcpp::NumericVector x = df["x"];
    Rcpp::NumericVector y = df["y"];
    Rcpp::NumericVector xx = x[ indices ];
    Rcpp::NumericVector yy = y[ indices ];

    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";
    Rcpp::String line_id = polygon_id;

    std::string xyzm = "XY";
    bool close = false;

    res[ i ] = sfheaders::sfg::sfg_multipolygon( df_tri, geometry_cols, polygon_id, line_id, xyzm, close );

  }

  sfheaders::sfc::attach_sfc_attributes( res, attributes );

  return res;

}
library(sf)
library(sfheaders)

nc <- sf::st_read( system.file("./shape/nc.shp", package = "sf"))
sf_poly <- sfheaders::sf_cast( nc, "POLYGON" )

res <- decido:::earcut_sfc( sf_poly$geometry )
sf_tri <- sf::st_as_sf( res )
sf_tri <- sfheaders::sf_cast( sf_tri, "POLYGON" )
sf_tri$id <- sample( 1:nrow( sf_tri), replace = FALSE )

library(mapdeck)

set_token( secret::get_secret( "mapbox" ) )

mapdeck() %>%
  add_polygon(
    data = sf_tri
    , fill_colour = "id"
  )

Screen Shot 2020-05-03 at 7 17 47 pm

@mdsumner
Copy link
Member Author

mdsumner commented May 3, 2020

I'm not sure, but I think it doesn't belong here. But then, it's such a handy facility to have ... probably there's no harm in it and maybe we find a different arrangement later. I'll think.

In sfdct I return a GEOMETRYCOLLECTION, so the result is one-to-one with the input sfdf. There is a TIN type (it still requires 4 vertices) and that might be worthwhile? (but there's no support in sf, no plot, no cast - easy to fix the plot problem with your cast api I suppose )

https://github.com/hypertidy/silicate/blob/master/R/TIN.R

There's also TRIANGLE but I think TIN is more appropriate: fwiw:

sf::st_as_sfc("TRIANGLE ((0.1 0.1, 0.1 1.1, 1.1 1.1, 0.1 0.1))")

@dcooley
Copy link
Contributor

dcooley commented May 3, 2020

Laying my cards on the table, my sole motivation for making triangles directly in C++ is to by-pass Deck.gl's polygon-tesselator, which itself uses the javascript version of earcut. If I can get all the triangles made, and interleave the coordinates, then I can bypass all deck.gl's worker functions and go straight from C++ to the WebGL buffers for drawing (with very minimal javascript in between).

So in summary, I don't mind if this functionality doesn't have a home here, but this is an example of how it would work nonetheless.

@mdsumner
Copy link
Member Author

mdsumner commented May 3, 2020

Cool, I think it doesn't go here then. With the right tools as you've shown it's easy to arrange the parts, and we're working out the landscape of what we end up with.

@mdsumner mdsumner removed the help wanted Extra attention is needed label May 5, 2020
@dcooley
Copy link
Contributor

dcooley commented May 6, 2020

I think I mentioned I was experimenting returning the triangle coordinates directly (rather than the indices), and I've had some success here

I'm not sure where the correct home for it is though. It would make logical sense to move it here to decido, but, I've edited the earcut.hpp code just enough to mean you (the package) would need two versions of earcut.hpp. Which may not be ideal.

@mdsumner
Copy link
Member Author

mdsumner commented May 6, 2020

Appreciate the update, I'll put some thought into it in the next few days. I don't really mind the problem of duplicated source code, as long as there's a clear pathway. Your forking of earcut.hpp reminded me that decido needs the updated version, and some process to trigger updates in future would be good too. You must have encountered this with other libs?

@dcooley
Copy link
Contributor

dcooley commented May 6, 2020

You must have encountered this with other libs?

yeah, but I have never got a "process" sorted out, other than re-basing my fork, copy & pasting the source code...

@mdsumner
Copy link
Member Author

mdsumner commented May 6, 2020

heh

@tim-salabim
Copy link

yeah, but I have never got a "process" sorted out, other than re-basing my fork, copy & pasting the source code...

Same here... Though I've seen that "dependabot" (or whatever it's called) more and more in other repos recently.

@mdsumner
Copy link
Member Author

My process now, having update the source for the first time since original release:

https://github.com/hypertidy/decido/tree/master/inst/earcut_source

sf will come from elsewhere ;)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants