[R] Get contour of a map

Pierre Bruyer pbruyer at agaetis.fr
Tue May 24 14:48:06 CEST 2011


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:

- a scale of colors for the rasters (and the colors).
- for aesthetic, the contours of each scale must be smooth

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
> 



More information about the R-help mailing list