r-lib / isoband

isoband: An R package to generate contour lines and polygons.
https://isoband.r-lib.org
Other
131 stars 14 forks source link

helper for degenerate matrix #5

Open mdsumner opened 5 years ago

mdsumner commented 5 years ago

I was caught out by my own assumptions in a wrapper I made because this gives a "degenerate", 1-column matrix.

library(raster)
as.matrix(brick(raster(volcano)))  

What I saw was this error, the dimensions of what I was giving isobands was wrong, so this message is accurate:

Error in isobands_impl(x, y, z, levels_low, levels_high) : Not a matrix.

But, when it's a single row or column should isoband catch this case [1, ] or [ ,1] and report back immediately?

Currently (81b278c):

library(isoband)
library(isoband)
x <- 1:10
y <- 1:5
m <- t(volcano[1:10, 1:5])
#isobands(x, y, m, 101, 106)
i <- 2
## both cases return empty results
isobands(x[i], y, m[, i, drop = FALSE], 101, 106)
isobands(x, y[i], m[i, , drop = FALSE], 101, 106)

Could this column or row matrix case ever return something useable?

(It's otherwise working fabbo on lots of real-world data and I can't fault it!).

clauswilke commented 5 years ago

The "Not a matrix" error is generated by Rcpp, I believe. I am not explicitly checking for non-matrix input, but the C++ function takes a matrix as an argument and Rcpp seems to make sure that it can actually provide the function with a matrix.

For the degenerate cases (single row/single column), I feel the current output is correct. It is the same as you would get if you had provided a non-degenerate matrix with NA values or used levels that are outside the data range. The advantage of doing it this way (rather than erroring out or returning NULL) is that you get the same defined behavior every time, regardless of the input. Consequently, you don't later have to worry whether specific levels exist or not in the returned object. They always exist, but they may be empty.

clauswilke commented 5 years ago

(Btw., I believe @tylermorganwall is still running all the images for his rayshaderbot through isoband, and things seem to be going fine, so that's a pretty solid real-world test. The code seems to work.)

tylermorganwall commented 5 years ago

I did run into errors on some images and had to wrap my code with a tryCatch call to handle those images outside of sf (doesn't sound like the same error, but it's an error nonetheless). Here's the example code:

library(magick)
library(magrittr)
library(ggplot2)
library(sf)
#> Linking to GEOS 3.6.1, GDAL 2.1.3, proj.4 4.9.3
library(isoband)
library(colorspace)

twitter_profile_url = "https://pbs.twimg.com/profile_images/937671987551522816/Xxn5L5Z8_400x400.jpg"

sf_from_image <- function(image) {
  image_gray <- image %>% image_quantize(colorspace = "gray")
  image_raster <- as.raster(image_gray)
  d <- dim(image_raster)
  m <- matrix(c((255-col2rgb(image_raster)[1,])), nrow = d[1], ncol = d[2], byrow = TRUE)
  b <- isobands(1:d[2], d[1]:1, m, 20*(-1:13), 20*(0:14))
  bands <- iso_to_sfg(b)
  data <- st_sf(
    level = letters[1:length(bands)],
    geometry = st_sfc(bands)
  )
}

img = twitter_profile_url %>%
  image_read()
img_sf = sf_from_image(img)
sf_gg = ggplot(img_sf) +
  geom_sf(aes(fill=level,color=level)) +
  coord_sf(expand = FALSE) +
  scale_fill_discrete_sequential(aesthetics = c("color","fill"), palette = "Viridis", guide = "none", rev = TRUE) +
  theme_void()
sf_gg

#> Error in CPL_geos_is_empty(st_geometry(x)): Evaluation error: IllegalArgumentException: Invalid number of points in LinearRing found 3 - must be 0 or >= 4.

Here are some examples that resulted in this error:

https://twitter.com/rayshaderbot/status/1095967829361393664 https://twitter.com/rayshaderbot/status/1093191703316627457 https://twitter.com/rayshaderbot/status/1093074816272396288

clauswilke commented 5 years ago

That seems to be related to this problem: https://twitter.com/mdsumner/status/1097817156123672576

A ring with 3 points makes logical sense (it is a line that goes forth and back), but if the specification doesn't allow it then we should just delete those cases in the conversion to sfg.

@mdsumner Can you think of simple input (say, a 5x5 matrix or so) that creates such a degenerate result, so we can add it to the test suite?

@tylermorganwall Was it always this one error or do you see any others?

tylermorganwall commented 5 years ago

It was always that error, although I didn’t do a comprehensive review of all the cases where it failed.

clauswilke commented 5 years ago

Over the next day or two I'll commit some code that skips invalid rings. Would you be willing/able to run the new code and see if any other errors crop up?

clauswilke commented 5 years ago

Ok, all of these seem to work now. I hope I didn't mess up anything else.

library(magick)
#> Linking to ImageMagick 6.9.9.39
#> Enabled features: cairo, fontconfig, freetype, lcms, pango, rsvg, webp
#> Disabled features: fftw, ghostscript, x11
library(magrittr)
library(ggplot2)
library(sf)
#> Linking to GEOS 3.6.1, GDAL 2.1.3, PROJ 4.9.3
library(isoband)
library(colorspace)

sf_from_image <- function(image) {
  image_gray <- image %>% image_quantize(colorspace = "gray")
  image_raster <- as.raster(image_gray)
  d <- dim(image_raster)
  m <- matrix(c((255-col2rgb(image_raster)[1,])), nrow = d[1], ncol = d[2], byrow = TRUE)
  b <- isobands(1:d[2], d[1]:1, m, 20*(-1:13), 20*(0:14))
  bands <- iso_to_sfg(b)
  data <- st_sf(
    level = letters[1:length(bands)],
    geometry = st_sfc(bands)
  )
}

draw_image <- function(twitter_profile_url) {
  img <- image_read(twitter_profile_url)
  img_sf <- sf_from_image(img)
  ggplot(img_sf) +
    geom_sf(aes(fill=level,color=level)) +
    coord_sf(expand = FALSE) +
    scale_fill_discrete_sequential(aesthetics = c("color","fill"), palette = "Viridis", guide = "none", rev = TRUE) +
    theme_void()
}

draw_image("https://pbs.twimg.com/profile_images/937671987551522816/Xxn5L5Z8_400x400.jpg")

draw_image("https://pbs.twimg.com/profile_images/1047564038174126080/gkpsqEAd_400x400.jpg")

draw_image("https://pbs.twimg.com/profile_images/1091232174827802624/FB9b1NRp_400x400.jpg")

Created on 2019-02-19 by the reprex package (v0.2.1)