[R-sig-Geo] gstat memory

[Edzer Pebesma]

Murray Richardson wrote:
> 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)))
I'd say it has to be a memory leak. Note that gstat was originally not
written as a library; most of it was written before 1995, when R had not
been released as open source. Finding memory leaks is a nice challenge,
but little rewarding now that RAM is so cheap.
>
> ... 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?
Routinely, and out of the box, for over 5 years. However, not under
windows. I even missed that windows already comes as a 64 bit port,
these days!

-- 
Edzer Pebesma
Institute for Geoinformatics (ifgi), University of Münster
Weseler Straße 253, 48151 Münster, Germany. Phone: +49 251
8333081, Fax: +49 251 8339763 http://ifgi.uni-muenster.de/
http://www.springer.com/978-0-387-78170-9 e.pebesma at wwu.de