윤곽선과 다각형을 교차 R

Question

나는 며칠 동안 등고선을 만들고 같은 파일에 윤곽선과 등고선을 그려 보려고 노력했다. 이제 같은 플롯에서 윤곽선과 shapefile을 만들 수 있습니다. 셰이프 파일로 컨투어를 클립하고 싶습니다. 셰이프 파일 만 보여줍니다. https://www.dropbox.com/sh/ztvmibsslr9ocmc/YOtiwB8p9p

스크립트 파일 image.scale.R는 다음 위치에서 찾을 수 있습니다 "https://www.dropbox.com/s/2f5s7cc02fpozk7/image.scale.R"

:

데이터 temp.csv

는 https://www.dropbox.com/s/mg2bo4rcr6n3dks/temp.csv Shape 파일은 다음 위치에서 찾을 수 있습니다이 링크를 찾을 수 있습니다 다음과 같이 내가 지금까지 사용하고

코드는 다음과 같습니다

## 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

다음과 같이 0

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) 

위의 코드의 출력은 다음과 같습니다

이

지금, 나는 다각형 Shape 파일에서 색의 그라데이션과 윤곽을 제거 만 교차 부분을 남겨 얻을 싶어요.

도움이 매우 감사합니다.

Research : Stack exchange Gis에서이 링크 https://gis.stackexchange.com/questions/25112/clip-depth-contour-with-spatial-polygon을 발견했으며이 방법을 따르려고했는데 윤곽을 생성하는 동안 항상 오류가 발생했습니다.

https://stat.ethz.ch/pipermail/r-sig-geo/2009-May/005793.html에 다른 유사한 스레드가 있습니다. 하지만 내 데이터 세트에서 작동하지 못했습니다.

이 시점에서 저를 도와 주신 Marc님께 감사 드리고 싶습니다.

감사합니다.

Answer 1

사실, @baptiste는 솔루션에 대한 강한 힌트 인 recent paper by Paul Murrell을 제공했습니다. Paul은 his personal website에서 얻을 수있는 전체 종이 원고 코드를 제공 해준 것에 관대했습니다. 측면 주제에서 Paul은이 논문을 통해 재현 가능한 연구의 아름다운 예를 보여줍니다. 일반적으로, 나는 다음과 같은 접근 방식에 나섭니다 : 추출 위도와 경도가 Shape 파일에서 좌표

• (이 작업을 수행하는 함수가 폴 Hiemstra에 의해, here이다), 코드와
• 플롯 모든 것을

• 및 polypath을 사용하여 추출 된 좌표를 기준으로 shapefile에 의해 정의 된 경계 외부의 모든 것을 제거하십시오.

  #function to extract coordinates from shapefile (by Paul Hiemstra) 
  allcoordinates_lapply = function(x) { 
      polys = x@polygons 
      return(do.call("rbind", lapply(polys, function(pp) { 
        do.call("rbind", lapply(pp@Polygons, coordinates)) 
        }))) 
  } 
  q = allcoordinates_lapply(shapefile) 
  
  #extract subset of coordinates, otherwise strange line connections occur... 
  lat = q[110:600,1] 
  long = q[110:600,2] 
  
  #define ranges for polypath 
  xrange <- range(lat, na.rm=TRUE) 
  yrange <- range(long, na.rm=TRUE) 
  xbox <- xrange + c(-20000, 20000) 
  ybox <- yrange + c(-20000, 20000) 
  
  #plot your stuff 
  plot(shapefile, lwd=2) 
  image(fld, axes=F, add=T) 
  contour(fld, add=T) 
  #and here is the magic 
  polypath(c(lat, NA, c(xbox, rev(xbox))), 
       c(long, NA, rep(ybox, each=2)), 
       col="white", rule="evenodd")