[R-sig-Geo] convert x,y,z data into a grid

[Meesters, Erik]

Michal,
thanks for the suggestions. If I try the wireframe I get a warning
message

"Warning message:
In nx * ny : NAs produced by integer overflow"

I haven't been able to figure out why that is, but could have something
to do with the 2.5mil * 2.5mil number of rows....

My data are not in a grid format, so the second option will not work. 
The data are simply UTM positions and depth values from a cruise
covering approximately an area of 10 square kilometers. In Grass there
is a function that reads in the dataframe and converts the data to a
grid which I think is easier to work with in many of the different
plotting routines.

Hope anyone can help.
Erik



-----Original Message-----
From: Michal Gallay [mailto:mgallay01 at qub.ac.uk] 
Sent: Wednesday, May 07, 2008 11:48
To: Meesters, Erik
Subject: Re: [R-sig-Geo] convert x,y,z data into a grid

On May 7 2008, Meesters, Erik wrote:

> Dear list members,
> I'm looking for a function in R like the function r.in.xyz in GRASS
(vs.
> 6.2). I have x, y, z data frame (UTM) with z equals depth and would 
> like to make a 3dperspective plot with rgl. My dataframe has 2.5 
> milion rows, so I think it's necessary to convert the data to a grid 
> with something like mean depth per grid cell. I know I've seen an 
> answer to a similar question somewhere on one of the lists, but I just

> don't seem to be able to find it again. Thanks for any help.
> Erik Meesters
>  

Hello Erik,

you could try something like the following code with lattice package.
You don't need to convert the dat to grid for this.

require(lattice)

mydata <- read.table(...)

wireframe(z ~ x * y, data = mydata,
          scales = list(arrows = FALSE),
          drape = TRUE, colorkey = TRUE,
          screen = list(z = 30, x = -60))

or with simple graphics function 'persp' but for this you need to
convert the data to matrix/grid. The following code presumes your x,y,z
coordinates are on a grid. It will not work if they are not or at least
I am not aware of any other approach.

require(sp)

coordinates(mydata) <- ~ x+y
gridded(mydata) <- TRUE
m <- as.matrix(mydata)
exag <- 3
spacing <- 5
theta_val <- 120
phi_val <- 35
ltheta_val <- 20
lphi_val <- 10
z <- exag * m      # Exaggerate the relief
x <- spacing * (1:dim(m)[1])
y <- spacing * (1:dim(m)[2])

par(c(0.5,0.5,0.5,0.5))
persp(x, y, z, theta = theta_val, phi = phi_val, col = "gray", scale =
FALSE, ltheta = ltheta_val, lphi=lphi_val, shade = 0.75, border = NA,
box = TRUE, r=600)

Both of them worked for my 2mil points data.

I hope this helps.

Michal


> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> 

--
Michal Gallay

Postgraduate Research Student
School of Geography, Archaeology and Palaeoecology Queen's University
Belfast BT7 1NN Northern Ireland

Tel: +44(0)2890 273929
Fax: +44(0)2890 973212
email: mgallay01 at qub.ac.uk
www: www.qub.ac.uk/geog