In R, it seems none of the currently available options (e.g. image) allow for fast real-time display of 2D raster graphics. I was interested for example to make a real-time interactive Mandelbrot viewer, where I was hoping to display 1920x1080 raw hex color matrices as raster images to achieve ca. 5-10 fps (calculating the Mandelbrot images themselves now achieves ca. 20-30 fps at moderate zooms). Or even better I would like to be able to pass an intensity pixel map to the graphics card & have the intensity to colour mapping also done using OpenGL acceleration. Using image() with option useRaster=TRUE, plot.raster or even grid.raster() doesn't cut it as displaying the raster image takes ca. 1/4 of a second (more than actually calculating the frame), so I am on the lookout for a more performant option, ideally using SDL or OpenGL acceleration.
I noticed that one should be able to call SDL, SDL_image and GL/OpenGL functions from R using the rdyncall foreign function interface package, which should have much better performance.
Although this package is archived on CRAN, an updated version is available here (and we are trying to get it back on CRAN)
library(remotes)
remotes::install_github("hongyuanjia/rdyncall")
library(rdyncall)
The SDL, SDL_image and SDL_mixer DLLs (version 1.2) (on Windows) have to be installed first from https://libsdl.org/release/, https://www.libsdl.org/projects/SDL_image/release/ and https://www.libsdl.org/projects/SDL_mixer/release/)(the 64 bit DLLs are to be put underR-4.2.1/bin/x64`) or using
source("https://raw.githubusercontent.com/Jean-Romain/lidRviewer/master/sdl.R")
On Ubuntu they can be installed using
sudo apt-get install libsdl1.2-dev libsdl-image1.2-dev libsdl-mixer1.2
Some demos of how to make SDL & OpenGL calls are available at https://dyncall.org/demos/soulsalicious/index.html (1980s computer-game style starfield, with music included).
Am I correct that with this package one should be able to display a 2D image raster using SDL & opengl acceleration? If so, has anyone any thoughts how to do this (I'm asking because I have no experience in using either SDL or OpenGL)?
To open a 1920 x 1080 SDL window I think I have to use
(gathered from some OpenGL examples and windowed.R script in https://dyncall.org/demos/soulsalicious/soulsalicious.tar.gz, fullscreen also possible, see fullscreen.R)
init <- function()
{
require(rdyncall)
dynport(SDL)
SDL_Init(SDL_INIT_VIDEO)
dynport(GL)
dynport(GLU)
dynport(SDL_image)
SDL_GL_SetAttribute(SDL_GL_RED_SIZE,8)
SDL_GL_SetAttribute(SDL_GL_GREEN_SIZE,8)
SDL_GL_SetAttribute(SDL_GL_BLUE_SIZE,8)
SDL_GL_SetAttribute(SDL_GL_DOUBLEBUFFER,1)
x_res <- 1920
y_res <- 1080
win <- SDL_SetVideoMode(x_res, y_res, 32,
SDL_HWSURFACE + SDL_OPENGL + SDL_DOUBLEBUF)
SDL_WM_SetCaption("SDL surface",NULL)
glEnable(GL_TEXTURE_2D)
# Set the projection matrix for the image
glMatrixMode(GL_PROJECTION)
glLoadIdentity()
x_min=1
x_max=x_res
y_min=1
y_max=y_res
glOrtho(x_min, x_max, y_min, y_max, -1, 1)
# Set the modelview matrix for the image
glMatrixMode(GL_MODELVIEW)
glLoadIdentity()
}
init()
I gather I should then set up some pixel transfer map using something like
glPixelMapfv(GL_PIXEL_MAP_I_TO_R, nb_colors, map_colors)
glPixelMapfv(GL_PIXEL_MAP_I_TO_G, nb_colors, map_colors)
glPixelMapfv(GL_PIXEL_MAP_I_TO_B, nb_colors, map_colors)
then create a buffer for the pixel data using another pixels <- glPixelMapfv call & draw the pixel data to the screen using glDrawPixels and swap the back and front buffers to display the image
using SDL_GL_SwapBuffers(win) and then wait for the user to close the window & then clean up using SDL_Quit() etc. Trouble is I have no OpenGL or SDL experience, so would anybody know how to carry out these last few steps? (I am using SDL version 1.2 here)
Some timings of non-OpenGL options which are too slow for my application:
# some example data & desired colour mapping of [0-1] ranged data matrix
library(RColorBrewer)
ncol=1080
cols=colorRampPalette(RColorBrewer::brewer.pal(11, "RdYlBu"))(ncol)
colfun=colorRamp(RColorBrewer::brewer.pal(11, "RdYlBu"))
col = rgb(colfun(seq(0,1, length.out = ncol)), max = 255)
mat=matrix(seq(1:1080)/1080,nrow=1920,ncol=1080,byrow=TRUE)
mat2rast = function(mat, col) {
idx = findInterval(mat, seq(0, 1, length.out = length(col)))
colors = col[idx]
rastmat = t(matrix(colors, ncol = ncol(mat), nrow = nrow(mat), byrow = TRUE))
class(rastmat) = "raster"
return(rastmat)
}
system.time(mat2rast(mat, col)) # 0.24s
# plot.raster method - one of the best?
par(mar=c(0, 0, 0, 0))
system.time(plot(mat2rast(mat, col), asp=NA)) # 0.26s
# grid graphics - tie with plot.raster?
library(grid)
system.time(grid.raster(mat2rast(mat, col),interpolate=FALSE)) # 0.28s
# base R image()
par(mar=c(0, 0, 0, 0))
system.time(image(mat,axes=FALSE,useRaster=TRUE,col=cols)) # 0.74s # note Y is flipped to compared to 2 options above - but not so important as I can fill matrix the way I want
# ggplot2 - just for the record...
df=expand.grid(y=1:1080,x=1:1920)
df$z=seq(1,1080)/1080
library(ggplot2)
system.time({q <- qplot(data=df,x=x,y=y,fill=z,geom="raster") +
scale_x_continuous(expand = c(0,0)) +
scale_y_continuous(expand = c(0,0)) +
scale_fill_gradientn(colours = cols) +
theme_void() + theme(legend.position="none"); print(q)}) # 11s
Five years after first having posted this question finally found the answer:
One efficient and convenient option to display a raster image in R is to use the nativeRaster format of the nara package - this is ca. 10x faster than using image. Ca. 10x faster still, one can use an SDL+OpenGl solution which can be called via the rdyncall foreign function interface, and this solution also allows displaying graphics fullscreen as opposed to just in a window.
Theoretically, this allows for 200+ fps raster display in 1920x1080 resolution. I am using an OpenGL Texture mapping approach - I haven't tried using an SDL_BlitSurface approach instead - this will give different timings still. And doing the intensity to colour mapping using SDL or OpenGL presumably would also be possible, though I couldn't immediately figure out how to do that (as I now do this beforehand, this would still be a big bottleneck). If anyone knows how to achieve that in pure OpenGL (maybe using glPixelMapfv as mentioned in the question above), please let me know. Using SFML, interfaced via Rcpp, instead of SDL should also be possible and should be faster than using SDL, but I haven't tried that either.
First we install the requird packages & calculate a nice example Mandelbrot fractal image :
# 0. load required packages ####
library(remotes)
# install SDL libraries (SDL, SDL_image & SDL_mixer, version 1.2)
# from https://libsdl.org/release/,
# https://www.libsdl.org/projects/SDL_image/release/ and
# https://www.libsdl.org/projects/SDL_mixer/release/)
# 64 bit DLLs are to be put under R-4.2.3/bin/x64
# you can use
source("https://raw.githubusercontent.com/Jean-Romain/lidRviewer/master/sdl.R")
# On Ubuntu install using sudo apt-get install libsdl1.2-dev libsdl-image1.2-dev libsdl-mixer1.2
remotes::install_github("hongyuanjia/rdyncall")
library(rdyncall)
remotes::install_github("coolbutuseless/nara")
library(nara)
remotes::install_github("tomwenseleers/mandelExplorer") # for nice example images
library(mandelExplorer)
library(microbenchmark)
# 1. create nice 640x480 image - here Mandelbrot fractal ####
xlims=c(-0.74877,-0.74872)
ylims=c(0.065053,0.065103)
x_res=640L
y_res=480L
nb_iter=as.double(nrofiterations(xlims))
m <- matrix(mandelRcpp2(as.double(xlims[[1]]), as.double(xlims[[2]]), # openmp+SIMD version
as.double(ylims[[1]]), as.double(ylims[[2]]),
x_res,
y_res,
nb_iter),
nrow=as.integer(x_res))
m[m==0] <- nb_iter
m <- equalizeman(m, nb_iter, rng = c(0, 0.95), levels = 1E4)^(1/8) # equalize colours & apply gamma correction
dim(m) # x_res x y_res, grayscale matrix normalized between 0 and 1
range(m) # 0 1
x11(type = 'dbcairo', antialias = 'none', width = 6*x_res/y_res, height = 6, bg = "black", canvas = "black") # Setup a fast graphics device that can render quickly
par(mar=c(0, 0, 0, 0))
microbenchmark(image(m, col=palettes[[2]], asp=y_res/x_res, axes=FALSE, useRaster=TRUE), unit='s')
# 0.075s = 13 fps, this includes colour mapping+raster display
x11(type = 'dbcairo', antialias = 'none', width = 6*x_res/y_res, height = 6, bg = "black", canvas = "black")
microbenchmark({natrast = matrix_to_nr(mat=matrix(m[,ncol(m):1], ncol = x_res, nrow = y_res, byrow = FALSE),
palette=palettes[[2]],
fill = "black")
plot(natrast);
dev.flush()}, unit='s')
# 0.007s=143 fps=11x faster than image, including the colour mapping + raster display
# timings just for display of nativeRaster
microbenchmark({ plot(natrast);
dev.flush() }, unit='s')
# 0.0043 s = 232 fps = 17x faster than image
dynport(SDL)
dynport(GL)
dynport(GLU)
# initialize SDL and create a window
if(SDL_Init(SDL_INIT_VIDEO) != 0) {
stop("Failed to initialize SDL")
}
# set OpenGL attributes
SDL_GL_SetAttribute(SDL_GL_RED_SIZE, 8)
SDL_GL_SetAttribute(SDL_GL_GREEN_SIZE, 8)
SDL_GL_SetAttribute(SDL_GL_BLUE_SIZE, 8)
SDL_GL_SetAttribute(SDL_GL_ALPHA_SIZE, 8)
SDL_GL_SetAttribute(SDL_GL_DOUBLEBUFFER, 1)
# set up screen: graphics window or run full screen
screen <- SDL_SetVideoMode(x_res, y_res, 0, SDL_OPENGL+SDL_DOUBLEBUF) # for windowed mode
# screen <- SDL_SetVideoMode(x_res, y_res, 0, SDL_OPENGL+SDL_DOUBLEBUF+SDL_FULLSCREEN) # for fullscreen mode
if(is.null(screen)) {
stop("Failed to set video mode")
}
# initialize projection matrix
glMatrixMode(GL_PROJECTION)
glLoadIdentity()
glOrtho(0, x_res, y_res, 0, -1, 1)
# initialize modelview matrix
glMatrixMode(GL_MODELVIEW)
glLoadIdentity()
# generate a texture
texId <- as.integer(0)
glGenTextures(1, texId)
glBindTexture(GL_TEXTURE_2D, texId)
# set texture parameters
glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_NEAREST)
glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_NEAREST)
# convert colour vector to matrix
colorsMatrix <- matrix(colors,
nrow = y_res,
ncol = x_res,
byrow = TRUE)
# create image matrix (created as nativeRaster to make sure it is a contiguous block of memory)
imageMatrix <- nara::raster_to_nr(colorsMatrix)
# timings just for nativeRaster display
microbenchmark({
# upload the image data to the texture
glTexImage2D(GL_TEXTURE_2D, 0, 4,
x_res, y_res, 0, GL_RGBA, GL_UNSIGNED_BYTE, imageMatrix)
# clear the window
glClear(GL_COLOR_BUFFER_BIT)
# enable textures
glEnable(GL_TEXTURE_2D)
# draw a quad with the texture
glBegin(GL_QUADS)
glTexCoord2f(0, 1); glVertex2f(0, 0)
glTexCoord2f(1, 1); glVertex2f(x_res, 0)
glTexCoord2f(1, 0); glVertex2f(x_res, y_res)
glTexCoord2f(0, 0); glVertex2f(0, y_res)
glEnd()
# disable textures
glDisable(GL_TEXTURE_2D)
# swap buffers to display the result
SDL_GL_SwapBuffers() },
unit='s')
# 0.0004s - 10x faster than nara solution above
# close window
SDL_Quit()
Output images are shown here.