Interpolate an irregular grid to a regular grid

I have an irregular grid that I need to convert to a regular grid in order to use the option image useRaster=TRUEfor graphic devices. I can do this on a small scale by converting an irregular grid to points, then interpolate the points with akima interp. However, it scales horribly with large sizes, so I'm looking for options.

Firstly, here is an example of a small scale (5x10), where only the x-dimension is irregular:

nx <- 5
ny <- 10
si <- list()  # irregular surface
si$x <- cumsum(runif(nx) * 10) + 100
si$y <- seq(20, 50, length.out=ny)
si$z <- matrix(rnorm(nx * ny), ncol=ny)
image(si)

And bilinear interpolated result:

sr_x <- seq(min(si$x), max(si$x), length.out=nx * 5)
sr_y <- si$y  # this dimension is already regular
require(akima)  # interpolate from points repeated off irregular grid
sr <- interp(rep(si$x, length(si$y)), rep(si$y, each=length(si$x)), si$z,
             xo=sr_x, yo=sr_y)
image(sr, useRaster=TRUE)

However, if a larger non-standard mesh is used (for example, nx <- 50; ny <- 100), the procedure is very slow. Is there a library or technique that will speed up the process?

, , .. ( ), 0,5 30 , 365 . , . , .

, :

nx <- 200
ny <- 10
si <- list()  # irregular surface
si$x <- cumsum(runif(nx, 0.5, 30) / 365)
si$y <- 1:ny
si$z <- matrix(rnorm(nx * ny), ncol=ny)
image(si)

:

dx <- 1/12  # 1 month spacing along x-axis
sr <- list()  # regular surface
sr$x <- seq(min(si$x), max(si$y), dx)  # equal-width breaks
nsrx <- length(sr$x)
sr$y <- si$y  # this dimension is already regular
sr$z <- matrix(nrow=length(sr$x), ncol=length(sr$y))
# Classify irregular dimension
si_xc <- cut(si$x, sr$x, include.lowest=TRUE, labels=FALSE)
# Aggregate means from irregular to regular dimension
for(xi in seq_len(nsrx))
    sr$z[xi,] <- apply(si$z[si_xc == xi, , drop=FALSE], 2, mean)
image(sr, zlim=range(si$z), useRaster=TRUE)

, , datsets 100s . , .