[R-sig-Geo] gstat memory

[Murray Richardson]

Hello R-SIG-GEO list,

I know this has come up before but I am having an ongoing memory problem 
with the gstat package (gstat out of dynamic memory) that I can't seem 
to solve.

I am using R to interpolate DEMs from LiDAR xyz point files and mosaic 
them together via RSAGA.  The script uses a loop to load in each xyz 
point file, and interpolate over a regular lattice of points using idw 
from gstat.  Although large, the computational requirements within each 
iteration of the loop should be well within my system's ability so it 
seems like it is a cumulative effect (note I can restart the process at 
the last loop that triggered the error and the iteration is 
successful).  I am removing temp objects and running gc() at the end of 
every loop.

It proceeds normally for about 5-10 iterations and physical memory use 
on my system (VISTA 64 bit, 12GB RAM) gradually increases over time 
until I get the "gstat out of memory" error. Here is the relevant 
portion of the loop:

...

   path<-"E:/LidarData/1_Ground_First_Return/UTM17/"
   filen<-paste(path, centre_tile$TILECODE, ".xyz", sep="")
   xyzi<-read.table(filen, sep=",", header=F)
   names(xyzi)<-c("x","y","z","i")

   for(j in 1:length(tile_group_names)){
       filen<-paste(path, tile_group_names[j], ".xyz", sep="")
       tmp_xyzi<-read.table(filen, sep=",", header=F)
       names(tmp_xyzi)<-c("x","y","z","i")
       xyzi<-rbind(xyzi,tmp_xyzi)

   }
       coordinates(xyzi)=~x+y
   # note grid_coords is just the regular lattice that is created from 
the current tile coordinates
   coordinates(grid_coords)=~x+y

   # do the interpolation
   interp<-idw(z~1, xyzi, grid_coords, nmax=4,maxdist=2, idp=1.0)

   names(interp)<-c("z","var")
   slot(interp,"data")<-data.frame(slot(interp,"data"))
   finalSPntsDF<-SpatialPointsDataFrame(interp, data.frame(interp$z),
               proj4string = CRS(as.character(NA)), match.ID = TRUE)


   finalSPDF<-SpatialPixelsDataFrame(finalSPntsDF, data.frame(interp$z), 
tolerance = sqrt(.Machine$double.eps),
     proj4string = CRS(as.character(NA)))

... 
On a related note - I have tried using the 64 bit version of Revolution 
but unfortunately the gstat package has not been ported.  Has anyone 
contemplated or begun a 64 bit port of this package?

Thanks

Murray