Crop dla SpatialPolygonsDataFrame

Question

Mam dwa pliki: SpatialPolygonsDataFrame DAT1, dat2

extent(dat1) 
class  : Extent 
xmin  : -180 
xmax  : 180 
ymin  : -90 
ymax  : 90 


extent(dat2) 
class  : Extent 
xmin  : -120.0014 
xmax  : -109.9997 
ymin  : 48.99944 
ymax  : 60 

chcę przyciąć DAT1 pliku przy użyciu zakresu dat2. Nie wiem, jak to zrobić. Po prostu obsługuję pliki rastrowe za pomocą funkcji "crop".

Kiedy używać tej funkcji dla moich bieżących danych, pojawia się następujący błąd:

> r1 <- crop(BiomassCarbon.shp,alberta.shp) 
Error in function (classes, fdef, mtable) : 

unable to find an inherited method for function ‘crop’ for signature"SpatialPolygonsDataFrame"’ 

Answer 1

sposób usprawnione dodano 2014-10-9:

raster::crop() mogą być stosowane w celu wycięcia Spatial* (jak również Raster*) obiekty.

Na przykład, oto jak można go używać w celu wycięcia SpatialPolygons* obiektu:

## Load raster package and an example SpatialPolygonsDataFrame 
library(raster) 
data("wrld_simpl", package="maptools") 

## Crop to the desired extent, then plot 
out <- crop(wrld_simpl, extent(130, 180, 40, 70)) 
plot(out, col="khaki", bg="azure2") 

---

Original (i nadal funkcjonalne) Odpowiedź:

W rgeos funkcja gIntersection() marki to całkiem proste.

Korzystanie ładne przykład mnel jako baza wypadowa:

library(maptools) 
library(raster) ## To convert an "Extent" object to a "SpatialPolygons" object. 
library(rgeos) 
data(wrld_simpl) 

## Create the clipping polygon 
CP <- as(extent(130, 180, 40, 70), "SpatialPolygons") 
proj4string(CP) <- CRS(proj4string(wrld_simpl)) 

## Clip the map 
out <- gIntersection(wrld_simpl, CP, byid=TRUE) 

## Plot the output 
plot(out, col="khaki", bg="azure2") 

Answer 2

zobaczyć upraw

corp (x, y, filename = "", snap = 'w pobliżu? ”Typ danych = NULL, ...)

x raster * przedmiotem

Y Zakres obiekt lub dowolny Obiekt, z którego obiekt stopniu możliwa jest wyodrębnione (szczegóły

Trzeba rasteryzację pierwszy SpatialPolygon korzystając rasterize funkcji z pakietu rastrowej

tworzę jakieś dane, aby pokazać, jak korzystać Rasteryzuj :

n <- 1000 
x <- runif(n) * 360 - 180 
y <- runif(n) * 180 - 90 
xy <- cbind(x, y) 
vals <- 1:n 
p1 <- data.frame(xy, name=vals) 
p2 <- data.frame(xy, name=vals) 
coordinates(p1) <- ~x+y 
coordinates(p2) <- ~x+y 

gdy próbuję:

crop(p1,p2) 
unable to find an inherited method for function ‘crop’ for signature ‘"SpatialPointsDataFrame"’ 

Teraz używając Rasteryzuj

r <- rasterize(p1, r, 'name', fun=min) 
crop(r,p2) 

class  : RasterLayer 
dimensions : 18, 36, 648 (nrow, ncol, ncell) 
resolution : 10, 10 (x, y) 
extent  : -180, 180, -90, 90 (xmin, xmax, ymin, ymax) 
coord. ref. : +proj=longlat +datum=WGS84 
data source : in memory 
names  : layer 
values  : 1, 997 (min, max) 

Answer 3

Nie można używać uprawy na obiektach sp wielokąta. Będziesz musiał utworzyć wielobok reprezentujący współrzędne bboxa dat2, a następnie użyć gIntersects w bibliotece rgeos.

Edytuj: Ten komentarz był związany z wersją dostępną w 2012 roku i nie jest to już prawdą.

Answer 4

Oto przykład jak to zrobić z rgeos używając mapy świata jako przykład

ten pochodzi z Rogerem Bivand na R-sig-Geo mailing list. Roger jest jednym z autorów pakietu sp.

Stosując mapę świata jako przykład

library(maptools) 

data(wrld_simpl) 

# interested in the arealy bounded by the following rectangle 
# rect(130, 40, 180, 70) 

library(rgeos) 
# create a polygon that defines the boundary 
bnds <- cbind(x=c(130, 130, 180, 180, 130), y=c(40, 70, 70, 40, 40)) 
# convert to a spatial polygons object with the same CRS 
SP <- SpatialPolygons(list(Polygons(list(Polygon(bnds)), "1")), 
proj4string=CRS(proj4string(wrld_simpl))) 
# find the intersection with the original SPDF 
gI <- gIntersects(wrld_simpl, SP, byid=TRUE) 
# create the new spatial polygons object. 
out <- vector(mode="list", length=length(which(gI))) 
ii <- 1 
for (i in seq(along=gI)) if (gI[i]) { 
    out[[ii]] <- gIntersection(wrld_simpl[i,], SP) 
    row.names(out[[ii]]) <- row.names(wrld_simpl)[i]; ii <- ii+1 
} 
# use rbind.SpatialPolygons method to combine into a new object. 
out1 <- do.call("rbind", out) 
# look here is Eastern Russia and a bit of Japan and China. 
plot(out1, col = "khaki", bg = "azure2")