[R-sig-Geo] convert x,y,z data into a grid
Meesters, Erik
Erik.Meesters at wur.nl
Wed May 7 14:23:09 CEST 2008
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
More information about the R-sig-Geo
mailing list