Wie erstelle ich ein räumliches Gitter aus einem Raster?

9

Ich benötige ein räumliches Raster als Master-Raster für verschiedene thematische Karten. Wie erstelle ich ein räumliches Gitter aus einem Raster, bei dem alle NA-Pixel verworfen werden?

Janbob Squarebrains
quelle
6
Wirf uns hier ein paar Fetzen. Ein kleiner Code zum Erstellen eines Rasters und was erwarten Sie davon?
Spacedman

Antworten:

14

Sie können alle Nicht-NA-Koordinaten von Zellen in einem Raster abrufen mit:

r = raster(matrix(runif(20),5,4))
r[r>.5]=NA

coordinates(r)[!is.na(values(r)),]
          x   y
 [1,] 0.375 0.7
 [2,] 0.125 0.5
 [3,] 0.375 0.5
 [4,] 0.625 0.5
 [5,] 0.875 0.5
 [6,] 0.125 0.3
 [7,] 0.375 0.3
 [8,] 0.625 0.3
 [9,] 0.375 0.1
[10,] 0.875 0.1

Das sind die Zellen, die nicht NA sind. Sie können diese dann wahrscheinlich SpatialPixels zuführen

SpatialPixels(SpatialPoints(coordinates(r)[!is.na(values(r)),]))
Object of class SpatialPixels
Grid topology:
  cellcentre.offset cellsize cells.dim
x             0.125     0.25         4
y             0.100     0.20         4
SpatialPoints:
          x   y
 [1,] 0.375 0.7
 [2,] 0.125 0.5
 [3,] 0.375 0.5
 [4,] 0.625 0.5
 [5,] 0.875 0.5
 [6,] 0.125 0.3
 [7,] 0.375 0.3
 [8,] 0.625 0.3
 [9,] 0.375 0.1
[10,] 0.875 0.1
Coordinate Reference System (CRS) arguments: NA 

Obwohl ich persönlich alles in einem Raster als Raster behalten würde.

Ich bin mir immer noch nicht ganz sicher, was Sie wollen - SpatialGridObjekte definieren vollständige rechteckige Gitter, daher macht eines ohne die NA-Pixel keinen Sinn.

Spacedman
quelle
danke Kumpel, ich bin mir auch nicht sicher, was ich will. Ich versuche, einen Workflow zu entwickeln, um aus verschiedenen Raumpunktdaten identische Rasterkarten (in Bezug auf Auflösung, Koordinaten, Crs, ...) zu erstellen. Rasterisierung ist aufgrund der räumlichen Punktstreuung keine Option. Müssen iwd oder ähnliches verwenden. Sie haben Recht, das Raster ist rechteckig - was ich als Master-Vorlage brauche, ist wahrscheinlich ein SpatialPointsDataFrame. was gerastert werden kann. S.
Janbob Squarebrains
Sie möchten wahrscheinlich ein Master-Raster (ha!) Und dann alle nachfolgenden Raster gleichzeitig erstellen.
Radar
2
Einverstanden. Wenn die Dinge gerastert sind, ist Raster fast immer die Antwort. Entweder das, oder Sie blättern in R über etwa 20 Seiten der angewandten Geodatenanalyse, um herauszufinden, wie SpatialPixel ausgeführt werden. Ich habe diese Wahl heute tatsächlich einem Doktoranden gegeben. Er wählte ... mit Bedacht!
Spacedman
20

Um einen RasterLayer in ein Spatial * -Objekt (Gitter oder Pixel) umzuwandeln, können Sie die Zwangsfunktion "as" verwenden.

library(raster)
r <- raster(matrix(runif(20),5,4))
r[r>.5] <- NA
g <- as(r, 'SpatialGridDataFrame')
p <- as(r, 'SpatialPixels')

plot(r)
points(p)
Robert Hijmans
quelle
7

Ihre beiden Anforderungen scheinen sich auf verschiedene Dinge zu beziehen:

1) Eine Art zuverlässige Raster-Raster-Vorlage.

2) Ein spärliches Gitter, in dem fehlende Zellen nicht explizit gespeichert werden.

sp :: GridTopology bietet die erste, es ist nur eine Beschreibung des Gitters basierend auf der unteren linken Zellkoordinate (cellcentre.offset), dem Zellabstand (Zellengröße) und den Abmessungen des Gitters (cells.dim).

Mit der sp :: SpatialPixelsDataFrame-Klasse können Sie spärliche Gitter speichern, aber allein speichert sie viel mehr als die "Vorlage" - sie speichert auch jede Koordinate explizit. Dies liegt daran, dass zwei Jobs ausgeführt werden. Zum einen können Sie die ursprünglichen Koordinaten beibehalten, die aus dem Raster stammen und möglicherweise leicht unregelmäßig sind. Zum anderen können Sie nur die Zellen speichern, die gültige Werte haben. (Diese beiden Ziele hätten wohl getrennt werden sollen, aber das ist eine andere Geschichte).

Ich glaube nicht, dass das Rasterpaket ein explizites Analogon zu GridTopology hat, aber Sie können die Komponenten beschaffen, um "Ihre eigenen zu rollen":

library(raster)
r1 <- raster(nrows=108, ncols=21, xmn=0, xmx=10)

## "cellsize"
res(r1)
## [1] 0.4761905 1.6666667

## extreme cell corners (just a different convention to sp's cellcentre)
extent(r1)
class       : Extent 
xmin        : 0 
xmax        : 10 
ymin        : -90 
ymax        : 90 

## we can also use bbox to get the same thing
bbox(r1)
min max
s1   0  10
s2 -90  90

## grid dimensions (including number of attributes/layers as 3rd "dim")
dim(r1)
## [1] 108  21   1

Wenn wir diese alle miteinander verbinden, können wir von Raster zu SP wechseln:

GridTopology(bbox(r1)[,1] + res(r1)/2, res(r1), dim(r1)[2:1])

(Beachten Sie, wie die Abmessungen umgekehrt werden müssen). Eine andere einfachere Möglichkeit besteht darin, auf SpatialGrid zuzugreifen und die getGridTopology von sp zu verwenden. Dies ist jedoch teurer, da die Koordinaten auf dem Weg generiert werden:

getGridTopology(as(r1, "SpatialGrid"))

Diese drei Teile der Raster- "Topologie" sind nicht alle erforderlich, da einige redundant sind, aber es gibt ansonsten keine Möglichkeit, sie alle als ein Objekt zu erfassen - außer, dass das oben erstellte Raster "leer" ist und dies tun kann Arbeit, die GridTopology für sp. Ich bin mir nicht sicher, wie "leer" es ist, aber es speichert den "Daten" -Slot sicherlich nicht explizit und ist kleiner als es mit Werten darin wäre. Das Rasterpaket bemüht sich im Allgemeinen, die Speichernutzung auf ein Minimum zu beschränken, sodass Sie sich möglicherweise keine Sorgen machen müssen, dass Sie wirklich "spärlich" sind.

Das könnte helfen, ein bisschen mehr zu erklären. Ich weiß, dass ich Spacedmans Antwort überlappe, aber es ist immer noch nicht klar, was Sie in der Frage genau meinen.

  • (Wahrscheinlich, da Sie spärliche Zellen speichern können, indem Sie nur ihren Index anstatt expliziter Koordinaten speichern, und ursprünglich "leicht unregelmäßige" Zellkoordinaten nur als Attribute gespeichert werden können, wenn Sie sie wirklich wollten. Weder sp noch Raster behandeln überhaupt unregelmäßige Raster. Nicht einmal der einfache geradlinige Fall - dies entspricht den meisten GIS-Tools, ist jedoch bedauerlich, da es in Formaten wie NetCDF durchaus üblich ist und von Rs guter alter Grafik :: Bild-Funktion verarbeitet wird (allerdings nicht, wenn die aktuelle Option rasterImage verwendet wird). .)
mdsumner
quelle
Vielen Dank für die sehr lehrreichen Aussagen, die ich verdauen muss. Hoffentlich kann ich mein Problem lösen oder eine explizitere Frage stellen!
Janbob Squarebrains
Auch @Spacedman: Ich habe meinen Fragenhintergrund hier im Forum
Janbob Squarebrains