hypertidy / decido

Constrained polygon triangulation by 'ear cutting'
https://hypertidy.github.io/decido
Other
15 stars 3 forks source link

sf prototype #9

Closed mdsumner closed 4 years ago

mdsumner commented 6 years ago

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 commented 6 years ago

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

dcooley commented 4 years ago

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

Things to consider

// [[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 commented 4 years ago

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 commented 4 years ago

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 commented 4 years ago

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.

dcooley commented 4 years ago

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 commented 4 years ago

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 commented 4 years ago

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 commented 4 years ago

heh

tim-salabim commented 4 years ago

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 commented 4 years ago

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 ;)