Re: [R] Get contour of a map

From: Pierre Bruyer <pbruyer_at_agaetis.fr>
Date: Tue, 24 May 2011 14:48:06 +0200

Yes sure! I understand english very well but I am not very good when I wrote, so I try to don't write too nonsense. I have post a precedent topic about that but no one have find a response for the ensemble of the problem.

So, I would like to create a raster for a map in R. I have a dataset with coordinates xy and for each points of a grid, and a z concentration in each points. I would obtain the same result than the software "Surfer 9". I have two constraints:

Now an exemple with imaginary data (My datasets have 10000 points):

x = c(0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0)
y = c(0,0,0,0,0,0,0,0,0,0,0,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.2,0.2,0.2,0.2,0.2,0.2,0.2,0.2,0.2,0.2,0.2,0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.4,0.4,0.4,0.4,0.4,0.4,0.4,0.4,0.4,0.4,0.4,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.6,0.6,0.6,0.6,0.6,0.6,0.6,0.6,0.6,0.6,0.6,0.7,0.7,0.7,0.7,0.7,0.7,0.7,0.7,0.7,0.7,0.7,0.8,0.8,0.8,0.8,0.8,0.8,0.8,0.8,0.8,0.8,0.8,0.9,0.9,0.9,0.9,0.9,0.9,0.9,0.9,0.9,0.9,0.9,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0)
z = c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,10,5,6,0,0,0,0,0,0,0,0,0,3,6,9,15,14,12,0,0,0,0,0,0,0,0,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0)

essai <- data.frame(x =x, y = y, z = z)

level3 <- matrix(0,6,4)

level3[1,1] <- 1
level3[2,1] <- 5
level3[3,1] <- 10
level3[4,1] <- 15
level3[5,1] <- 25
level3[6,1] <- 30

level3[1,2] <- 0
level3[2,2] <- 0
level3[3,2] <- 0
level3[4,2] <- 0

level3[5,2] <- 255
level3[6,2] <- 255
level3[1,3] <- 0
level3[2,3] <- 0
level3[3,3] <- 255
level3[4,3] <- 255
level3[5,3] <- 0
level3[6,3] <- 0

level3[1,4] <- 0
level3[2,4] <- 255
level3[3,4] <- 0
level3[4,4] <- 255

level3[5,4] <- 0
level3[6,4] <- 255

## grd is the xyz dataframe, level is the matrix level3,and map_output is the path where ##you want save the image.
create_map <- function(grd, level ,map_output, format = c("jpeg"), width_map = 150, height_map = 150,...) {

        grd2 <- matrix(grd[,3], nrow = sqrt(length(grd[,3])), ncol = sqrt(length(grd[,3])), byrow = FALSE)

	##creation of breaks for colors
	i<-1	
	paliers <- c(-1.0E300)
	while(i<=length(level[,1]))
	{
		paliers <- c(paliers,level[i,1])
		i <- i+1

}
paliers <- c(paliers, 1.0E300) ##scale color creation i <- 1 colgraph <- c(rgb(255,255,255, maxColorValue = 255)) while(i<=length(level[,2])) { colgraph <- c(colgraph, rgb(level[i,2],level[i,3],level[i,4], maxColorValue = 255)) i <- i +1
}

### now many possibilities, this is the possibility with get contour of map         

	grd2 <- rbind(-1E100, grd2, -1E100)	
	grd2 <- cbind(-1E100, grd2, -1E100)

	x <- seq(from=min(grd[,1]), to=max(grd[,1]), length = sqrt(length(grd[,1])))
	y <- seq(from=min(grd[,2]), to=max(grd[,2]), length = sqrt(length(grd[,2])))
	
	x <- c(min(grd[,1])-0.01, x , max(grd[,1])+0.01)
	y <- c(min(grd[,2])-0.01, y , max(grd[,2])+0.01)

	contours <- contourLines(x, y,grd2,levels=paliers)

	
	##user can choose the output format (default is jpeg)
	switch(format,  
		png = png(map_output, width = width_map, height = height_map) ,
		jpeg = jpeg(map_output, width = width_map, height = height_map, quality = 100),
		bmp = bmp(map_output, width = width_map, height = height_map),
		tiff = tiff(map_output, width = width_map, height = height_map),
		jpeg(map_output, width = width_map, height = height_map))


	## drawing map
		par(mar=c(0,0,0,0),xaxs="i", yaxs="i")

	plot(0,col="white",main="polygon()", asp = 1, axes = FALSE, ann = FALSE,xlim=c(min(grd[,1	]),max(grd[,1])), ylim = c(min(grd[,2]),max(grd[,2])),type = "n")

	for (i in seq_along(contours)) {
		 x <- contours[[i]]$x
		 y <- contours[[i]]$y
		 c <- contours[[i]]$level
		 j <- 1
		 tmp <- 0 
		 while(j < length(level[,1])+1 && tmp == 0){
			if( c == paliers[j]){
		 		tmp <- j
		 	}
		 	j <- j+1 
		}	
		a <- spline( seq_along(x), x)$y
		b <- spline( seq_along(y), y)$y
		polygon( a, b ,col = colgraph[tmp], border = NA)

}
dev.off()

}

Here I use the solution of Duncan Murdoch, it is almost that I want, now I try to resolve the default of his methods, e.g. I don't want smooth the contour ( the border) of my raster.

Thank for your help :-))
(Sorry for the spam, I have forgotten to send the mail to r-help)

Le 24 mai 2011 à 12:23, David Winsemius a écrit :

> 
> On May 24, 2011, at 5:40 AM, Pierre Bruyer wrote:
> 
>> I have seen this package and especially the 'perimeter' function, and it has not a 'levels' parameters like a lot of other functions, and this is a big problem for my project.
> 
> Without a reproducible example it is difficult to know what you expect from a `levels` parameter, but I can get particular contour levels to be drawn for crossed-cubic regression splines. It would be more productive if you offered a minimal example to be worked on.
> 
> -- 
> David.
>> Otherwise I find this package have a better architecture than the other spatial package ^^
>> 
>> 
>> Le 24 mai 2011 à 10:52, David Winsemius a écrit :
>> 
>>> 
>>> On May 23, 2011, at 5:55 AM, Pierre Bruyer wrote:
>>> 
>>>> Hello everybody,
>>>> 
>>>> I search a function which returns the contour of map with levels like contourLines, but I would like this function return the border of the map too, because the function contourLines cannot consider the corner of the map and it is not adapted to fill polygon after that.
>>> 
>>> Frank Harrell has a `perimeter` function in his `rms` package which he uses to establish (and optionally draw)  boundaries around 2D regions with sufficient data to yield meaningful estimates. The plotting is handed off to lattice functions. His default plotting function for 2D plotting is contourplot. It didn't take that long to make the transition to rms+lattice and I have been very pleased with the integration of the two.
>>> 
>>> There is also a chull  (convex hull) function although at the moment I do not remember which package it resides in.
>>> 
>>> 
>>>> 
>>>> Thanks in advance
>>>> 
>>>> Pierre Bruyer
>>>> ______________________________________________
>>>> R-help_at_r-project.org mailing list
>>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>>>> and provide commented, minimal, self-contained, reproducible code.
>>> 
>>> David Winsemius, MD
>>> West Hartford, CT
>>> 
>> 
> 
> David Winsemius, MD
> West Hartford, CT
> 

______________________________________________
R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Tue 24 May 2011 - 13:10:48 GMT

This quarter's messages: by month, or sorted: [ by date ] [ by thread ] [ by subject ] [ by author ]

All messages

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Tue 24 May 2011 - 13:20:08 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive