leaflet icon indicating copy to clipboard operation
leaflet copied to clipboard

leaflet with both stars RGB and sf polygon in custom CRS (epsg:3031)

Open CCheCastaldo opened this issue 2 years ago • 4 comments

I am trying to use R leaflet to display a stars RGB with a second layer that is a vector-based (spatial polygon sf). It is important that both spatial layers remain in polar stereographic (EPSG:3031) instead of having them reprojected to the standard EPSG:3857 web mercator projection used by leaflet. Also, I don't need to have any background tiles. When I do this I either get output that does not show the stars object if I specify the crs as 3031 but displays the polygon correctly in 3031 or, if I disable the custom crs, displays the stars object (which weirdly does not look reprojected to 3857) and the polygon (which is clearly reprojected to 3857). Again, what I am trying to do is have both layers be in 3031. Is this possible?

Here is a reproducible example:

library(sf)
library(terra)
library(stars)
library(leaflet)
library(leafem)

# create spatraster and then stars obj in 3031
coords <- c(xmin = -2399025, xmax = -2397255, ymin = 1702935, ymax = 1705005)
box <- terra::ext(coords)
r <- terra::rast(box, crs = "epsg:3031", resolution = 30)
r$swir1 <- rnorm(4071, 3, 1)
r$nir <- rnorm(4071, 2, 1)
r$red <- rnorm(4071, 1, 1)
rgb_stars <- stars::st_as_stars(r)

# create sf polygon obj in 3031
easting <- c(-2397904, -2398452, -2397612, -2397904)
northing <- c(1703249, 1704174, 1703763, 1703249)
drone_area = sf::st_sfc(sf::st_polygon(list(cbind(easting, northing))), crs = 3031)

# define the 3031 crs
crs_3031 <- leaflet::leafletCRS(
  crsClass = 'L.Proj.CRS',
  code = 'EPSG:3031',
  proj4def = '+proj=stere +lat_0=-90 +lat_ts=-71 +lon_0=0 +k=1 
  +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs',
  resolutions = 2^(16:1))

The leaflet generated here does not display the stars RGB at all but does display the polygon in 3031.

leaflet::leaflet(options = leaflet::leafletOptions(crs = crs_3031)) %>%
  addPolygons(data = drone_area %>% sf::st_transform(4326), color = "yellow") %>%
  addStarsRGB(rgb_stars, r = 3, g = 2, b = 1, group = "true color") %>%
  addLayersControl(overlayGroups = c("sites", "true color"))
Screen Shot 2022-05-21 at 9 50 18 PM

The leaflet generated here displays both the stars RGB and polygon but the polygon is, of course, transformed to 3857 while the stars RGB oddly looks like it is still in 3031.

leaflet::leaflet() %>%
  addPolygons(data = drone_area %>% sf::st_transform(4326), color = "yellow") %>%
  addStarsRGB(rgb_stars, r = 1, g = 2, b = 3, group = "true color") %>%
  addLayersControl(overlayGroups = c("sites", "true color"))
Screen Shot 2022-05-21 at 9 52 08 PM

Thanks for any help. I appreciate it.

CCheCastaldo avatar May 22 '22 01:05 CCheCastaldo

try {plainview} 🙏

mdsumner avatar Jul 09 '22 08:07 mdsumner

plainview doesn't do overlays... This is still an unsolved issue anywhere. I am hoping to get this done with georaster-for-leaflet.

tim-salabim avatar Jul 09 '22 08:07 tim-salabim

🙏

mdsumner avatar Jul 09 '22 08:07 mdsumner

plainview doesn't do overlays... This is still an unsolved issue anywhere. I am hoping to get this done with georaster-for-leaflet.

@tim-salabim many thanks for potentially addressing this in georaster-for-leaflet. I appreciate it. This issue has us dead in our tracks these past few months.

CCheCastaldo avatar Jul 29 '22 02:07 CCheCastaldo