Extrahieren Sie Raster aus Raster mit Polygon-Shapefile in R

Ich bin neu in R und benutze das Raster-Paket. Beim Extrahieren von Polygonen aus einer vorhandenen Rasterdatei ist ein Problem aufgetreten. Wenn ich benutze

extract(raster, poly_shape)

Funktion auf dem Raster erstellt es immer eine Liste mit den Daten. Ich möchte wirklich eine weitere Raster-Datei extrahieren, die ich wieder mit ArcGIS laden kann. Nachdem ich ein bisschen mehr gelesen habe, denke ich, dass die Zuschneidefunktion genau das ist, was ich wirklich brauche. Aber wenn ich versuche, diese Funktion zu nutzen

crop(raster, poly_shape)

Ich bekomme diesen Fehler:

Error in .local(x, y, ...) : extents do not overlap
In addition: Warning message:
In intersect(extent(x), extent(y)) : Objects do not overlap

Die Dateien raster und poly_shape sind für beide Funktionen gleich. Können Sie mir sagen, was hier falsch sein könnte? Ist es überhaupt richtig, dass die Zuschneidefunktion ein anderes Raster und keine Liste erstellt?

BEARBEITEN : Die Funktion Extent () funktioniert bei mir nicht. Ich bekomme immer noch den gleichen Fehler. Ich bin mir aber sicher, dass sich die beiden Datensätze überschneiden! Mit dem

extract(raster, poly_shape)

Ich bekomme die richtigen Daten daraus. Nur als Liste und nicht als Raster, wie ich es haben möchte. Ich habe gerade die Datasets in ArcGIS geladen, und sie passen sehr gut, sodass ich die Projektion nicht überprüft habe. Jetzt habe ich es versucht

projection(raster) # "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs"
projection(poly_shape) # "+proj=utm +zone=32 +ellps=GRS80 +units=m +no_defs"

und Sie können sehen, dass die Projektionen nicht passen. Die Extraktionsfunktion scheint in der Lage zu sein, die Dateien automatisch auf die richtige Weise zu transformieren. Ich weiß das, weil ich Folgendes getan habe:

• Ich habe den genauen Teil des Polygons, das ich in R extrahiert habe, auch in ArcGIS ausgeschnitten
• Ich habe die Summe aller Werte des extrahierten R-Polygons berechnet (Liste)
• Ich habe die Summe aller Rasterzellen berechnet, die ich in ArcGIS ausgeschnitten habe

Die 2 haben genau das gleiche Ergebnis, daher sollte die Schlussfolgerung lauten, dass die Extraktionsfunktion korrekt funktioniert hat. Jetzt habe ich zwei Möglichkeiten:

1. Ich brauche eine Möglichkeit, um einen Raster wieder aus der extrahierten Liste zu entfernen oder
2. Die beiden Datensätze (Raster + Poly_Shape) müssen dieselbe Projektion verwenden und die Zuschneidefunktion sollte funktionieren

Was würden Sie hier vorschlagen?

Die Extraktionsfunktion verhält sich genau so, wie sie sollte. Sie können erzwingen, dass die Zuschneidefunktion die Ausdehnung des Polygons verwendet, und dann das Objekt maskieren, um das genaue Raster zurückzugeben, das den Polygonbereich darstellt. Wenn Sie den Fehler weiterhin erhalten, bedeutet dies, dass sich Ihre Daten tatsächlich nicht überschneiden.

Bitte denken Sie daran, dass R keine "on the fly" -Projektion durchführt. Überprüfen Sie daher Ihre Projektionen. Mit der Funktion "extend ()" können Sie überprüfen, ob sich Ihre Extents überschneiden.

Hier ist ein Beispiel für das Beschneiden mithilfe der Polygonausdehnung und das anschließende Maskieren des resultierenden Rasters mithilfe des "gerasterten" Polygons.

# Add required packages
require(raster)
require(rgeos)
require(sp)

# Create some data using meuse 
data(meuse)
  coordinates(meuse) <- ~x+y
    proj4string(meuse) <- CRS("+init=epsg:28992")
data(meuse.grid)
  coordinates(meuse.grid) = ~x+y
    proj4string(meuse.grid) <- CRS("+init=epsg:28992")
      gridded(meuse.grid) = TRUE    
        r <- raster(meuse.grid) 
          r[] <- runif(ncell(r))

# Create a polygon
f <- gBuffer(meuse[10,], byid=FALSE, id=NULL, width=250, 
                         joinStyle="ROUND", quadsegs=10)   

# Plot full raster and polygon                       
plot(r)
  plot(f,add=T)

# Crop using extent, rasterize polygon and finally, create poly-raster
#          **** This is the code that you are after ****  
cr <- crop(r, extent(f), snap="out")                    
  fr <- rasterize(f, cr)   
    lr <- mask(x=cr, mask=fr)

# Plot results
plot(lr)
  plot(f,add=T)

Was ich eigentlich gesucht habe, war die mask()Funktion.

mask(raster, poly_shape)

funktioniert fehlerfrei und gibt das zurück, wonach ich gesucht habe.

Die Ausdehnung funktioniert einwandfrei ... Ich denke, die Xmin, Xmax, Ymin und Ymax Ihrer Ausdehnung unterscheiden sich von den X- und Y-Werten Ihres Rasters, dh sie sind entgegengesetzt gesetzt