R ggplot2 and Stars Package – Troubleshooting Slow Reprojection

ggplot2rstars

I'm finding very different plotting results with stars objects after I reproject them and use ggplot2. Consider the following:

library(stars)
library(ggplot2)

system.file("tif/L7_ETMs.tif", package = "stars") %>% read_stars() -> x

ggplot() + 
  geom_stars(data = x)

This works great, and is relatively quick. However, after I reproject

#add a transformation
x_proj <- st_transform(x, crs = 4326)

The next time I plot takes a looooong time. And for some plots (I was attempting to layer two rasters when I fell on this), it just fails. After running the above, try

ggplot() + 
  geom_stars(data = x_proj)

What is going on here? How can I restore the normal speed?

Best Answer

When you do that transformation...

> x_proj
stars object with 3 dimensions and 1 attribute
attribute(s):
             Min. 1st Qu. Median     Mean 3rd Qu. Max.
L7_ETMs.tif     1      54     69 68.91242      86  255
dimension(s):
     from  to offset delta refsys point                          values x/y
x       1 349     NA    NA WGS 84 FALSE [349x352] -34.9165,...,-34.8261 [x]
y       1 352     NA    NA WGS 84 FALSE  [349x352] -8.0408,...,-7.94995 [y]
band    1   6     NA    NA     NA    NA                            NULL    
curvilinear grid

Notice the last thing it says... "curvilinear grid". Its not a simple raster grid anymore, oh no. Each cell is now not rectangular. Its edges are not straight. To draw it, R has to transform every corner point separately, and work out the curve of each side of the cell, then draw all those 122,848 polygons.

If you warp it instead, you get back a regular grid on the new coordinate system, which plots quickly, because its a regular grid on the new coordinate system:

> x_warp = st_warp(x, crs = 4326)
> ggplot() +   geom_stars(data = x_warp)
> x_warp
stars object with 3 dimensions and 1 attribute
attribute(s):
             Min. 1st Qu. Median     Mean 3rd Qu. Max. NA's
L7_ETMs.tif     1      54     69 68.93072      86  255 8778
dimension(s):
     from  to   offset        delta refsys point values x/y
x       1 350 -34.9166  0.000259243 WGS 84    NA   NULL [x]
y       1 352 -7.94982 -0.000259243 WGS 84    NA   NULL [y]
band    1   6       NA           NA     NA    NA   NULL

That might be what you want, and you need to be aware of how warping may have modified your data (each cell could be a function of several cells, depending on how the warped CRS matches the original) and there could be smoothing or nearest-source-pixel mapping of the destination pixels.

If you do need the increased precision of the curvilinear transformed grid then I'm not sure there's much you can do to speed it up.

Related Question