Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

Hi I am wondering how can I create the contour lines from the discrete x, y and z dataset and then add the shapefile to cut the contour lines by polygon.

I want to clip the contour map by the shape file. I have created the contour map and added the shape file but I am not able to cut the shape file with polygon.

I used akima package to create the contour map. I have posted this question in stackoverflow R help but I didn't get any help so far and I desperately need help on this issue.

I will post my code and data links here :

The script file image.scale.R can be found on the following location "https://www.dropbox.com/s/2f5s7cc02fpozk7/image.scale.R "

The code that I have used so far is as follows:

## Required packages 
library(maptools) 
library(rgdal) 
library(sp) 
library(maptools) 
library(sm) 
require(akima) 
require(spplot) 
library(raster) 
library(rgeos)

## Set Working Directory 
setwd("C:\\Users\\jdbaba\\Documents\\R working folder\\shape")

## Read Data from a file 
age2100 <- read.table("temp.csv",header=TRUE,sep=",")

x <- age2100$x 
y <- age2100$y   
z <- age2100$z

####################################
##Load the shape file 
##################################### 
shapefile <- readShapePoly("Export_Output_4.shp")

fld <- interp(x,y,z)

par(mar=c(5,5,1,1)) filled.contour(fld)

###Import the image.scale 
source source("image.scale.R")

# http://menugget.blogspot.de/2011/08/adding-scale-to-image-plot.html

x11(width=8, height=7) 
layout(matrix(c(1,2), nrow=1, ncol=2), widths=c(6,1), height=6, respect=TRUE) 
layout.show(2)

par(mar=c(4,4,1,2)) 
image(fld,axes=T) 
contour(fld, add=TRUE)
#points(age2100$x,age2100$y, pch=".", cex=2,legend=F) 
plot(shapefile,add=T,lwd=2) 
box()
par(mar=c(4,0,1,4)) 
image.scale(fld, xlab="Eastings", ylab="Northings", xaxt="n", yaxt="n", horiz=FALSE)

axis(4)
mtext("Salinity", side=4, line=2.5)

The output of the above code is as follows:

![enter image description here][1]

Now, I want to get rid of the colored gradients and the contours from the polygon shapefile and only leave the intersection part.

Any help is highly appreciated.

Thanks.

share|improve this question
    
The data temp.csv can be found on this link dropbox.com/s/mg2bo4rcr6n3dks/temp.csv Shapefile can be found on the following location: dropbox.com/sh/ztvmibsslr9ocmc/YOtiwB8p9p –  Jdbaba Jan 18 '13 at 20:14
    
Why don't you just mask the interpolated data to the polygon before extracting contours from it? –  whuber Jan 18 '13 at 22:00
    
Thank you so much for your help. I don't know how to do it. Would you please give me suggestion how this can be performed ? I will try to follow your suggestion. Example of any setting would be helpful. –  Jdbaba Jan 18 '13 at 22:08

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.