Grid zuschneiden

Allgemeine Statistik mit R, die Test-Methode ist noch nicht bekannt, ich habe noch keinen Plan!

Moderatoren: EDi, jogo

Antworten
Everlast

Grid zuschneiden

Beitrag von Everlast »

Hallo (:

nachdem mir hier das letzte Mal bereits gut geholfen wurde bin ich in meinem Projekt inzwischen etwas weiter. Dabei ist eine neue Frage aufgekommen an der ich ein wenig fest hänge, war mir allerdings nicht sicher in welchem Teil des Forums ich es am besten poste. Zunächst einmal habe ich einen kleinen veränderten Beispieldatensatz. Es handelt sich dabei um Messwerte mit Koordinaten der jeweiligen Messstelle.

Code: Alles auswählen

dput(t1)
structure(list(Lfd = 1:20, MstNr = structure(c(16L, 17L, 15L, 
14L, 8L, 7L, 20L, 6L, 5L, 12L, 19L, 9L, 11L, 10L, 4L, 1L, 2L, 
18L, 13L, 3L), .Label = c("3369B", "3370", "3376", "3377", "3378", 
"3379", "3380", "3381", "3382", "3383", "3384", "3385", "3404", 
"3445", "3446", "3452", "3456", "636", "850", "851"), class = "factor"), 
    x = c(3469489.45369692, 3469083.73148078, 3469041.96280914, 
    3468635.21639453, 3468449.03013263, 3468450.69908635, 3469320.1627802, 
    3469424.6028951, 3469282.52506032, 3468388.84707303, 3469333.39947069, 
    3468243.43041729, 3469290.56665505, 3469588.87341771, 3469530.48687964, 
    3469971.09926508, 3469983.19236769, 3469736.45998308, 3470503.41533231, 
    3469088.22260067), y = c(5413581.63214282, 5413280.51059682, 
    5412268.69015278, 5412768.89058611, 5413137.30120601, 5412192.65287435, 
    5412100.59727255, 5414934.89892535, 5411931.3317301, 5413421.79058764, 
    5412878.37172573, 5413334.97727913, 5412852.64747422, 5413525.19319656, 
    5415130.55043546, 5412607.6942914, 5413845.45738916, 5412053.95384853, 
    5413575.24603885, 5415347.22405993), Messwert = c(0.333433333333333, 
    0.333433333333333, 0.334108333333333, 0.333433333333333, 
    0.356583333333333, 0.334401666666666, 0.351333333333333, 
    0.389209649122807, 1.15468717948718, 0.366882745098039, 0.443333333333333, 
    0.368706060606061, 0.35508696969697, 0.430333333333333, 0.497765151515151, 
    0.334573333333333, 0.366733833333333, 0.336677333333333, 
    0.392514583333333, 0.382595913978494)), .Names = c("Lfd", 
"MstNr", "x", "y", "Messwert"), class = "data.frame", row.names = c("1", 
"2", "3", "4", "5", "6", "7", "8", "9", "10", "11", "12", "13", 
"14", "15", "16", "17", "18", "19", "20"))
Was ich tun möchte ist diese auf ein Raster zu plotten und dieses Raster um die Punkte herum abzuschneiden. Den ersten Teil bekomme ich auch ganz gut hin denke ich:

Code: Alles auswählen

#Pakete 

install.packages("sp")
library(sp)

##Raumbezug herstellen
coordinates(t1) <- c("x", "y") 	# x/y in Koordinaten umwandeln
plot(t1)			# testweise plotten


###Raster dimensionieren

min_x = min(t1$x) - 100			# Min x + innerer Rand
min_y = min(t1$y) - 100			# Min y + innerer Rand
max_x = max(t1$x) 			# Max x			
max_y = max(t1$y) 			# Max y
x_length = max(t1$x - min_x) + 100	#  x-Länge + äußerer Rand
y_length = max(t1$y - min_y) + 100	#  y-Länge + äußerer Rand
cellsize = 100				# Zellgröße
ncol = round(x_length/cellsize,0)	# Zeilen
nrow = round(y_length/cellsize,0)	# Reihen

proj4string(t1) = CRS("+proj=utm +ellps=WGS84 +datum=WGS84")	# Projektion

#### Grid erstellen und Klasse zuweisen
grid = GridTopology(cellcentre.offset=c(min_x,min_y),cellsize=c(cellsize,cellsize),cells.dim=c(ncol,nrow))
grid = SpatialPixelsDataFrame(grid,
                              data=data.frame(id=1:prod(ncol,nrow)),
                              proj4string=CRS(proj4string(t1)))
                              
##### Punkte auf Grid plotten

plot(grid)
points(t1)
Leider hab ich keine wirklich Vorstellung wie ich die überschüssigen Teile des Grids nun abschneide. Ich hätte vermutet das man zunächst ein Polygon um die Punkte erstellt welches immer die äußersten Punkte verbindet. Anschließend könnte man das Polygon dazu benutzen die Fläche aus dem Raster zu clipen. Ich würde mich über jede Hilfe freuen (:

mit freundlichen Grüßen

Johann
Benutzeravatar
EDi
Beiträge: 1599
Registriert: Sa Okt 08, 2016 3:39 pm

Re: Grid zuschneiden

Beitrag von EDi »

Andere Idee:

Das raster direkt auf die richtige Größe (=bounding box + Rand) erzeugen:

Code: Alles auswählen

library(sp)
library(raster)
# make spatial
t1_sp <- t1 
coordinates(t1_sp) <- c("x", "y")
proj4string(t1_sp) = CRS("+proj=utm +ellps=WGS84 +datum=WGS84")	
plot(t1_sp)


# get bounding box of point & ad margin of 100
ext <- bbox(t1_sp)
ext[ , 1] <- ext[ , 1] - 100
ext[ , 2] <- ext[ , 2] + 100

# create raster from bounding box
r <- raster(extent(ext), 
            crs =  proj4string(t1_sp),
            resolution=100)
r
# add some random noise (just for color...)
r[] <- rnorm(length(r))

plot(r)
points(t1_sp)
Bitte immer ein reproduzierbares Minimalbeispiel angeben. Meinungen gehören mir und geben nicht die meines Brötchengebers wieder.

Dieser Beitrag ist lizensiert unter einer CC BY 4.0 Lizenz
Bild.
Everlast

Re: Grid zuschneiden

Beitrag von Everlast »

Hi EDI,

danke erstmal für den Vorschlag! Allerdings bekomme ich wenn ich das umsetze noch immer ein rechteckiges Grid. Ich stelle es mir aber eher so vor dass das Grid nur die Punkte umgibt, als würde man eine Linie um die Punkte herum zeichnen.

Vlt. kann ich es mit einem Beispiel etwas verdeutlichen, auch wenn hier teilweise ein etwas größerer Rand besteht und ich vermute dass die Kontur etwas anders zustande gekommen ist:

Code: Alles auswählen

library(gstat)

data(meuse.grid)
coordinates(meuse.grid) <- c("x", "y")
meuse.grid <- as(meuse.grid, "SpatialPixelsDataFrame")
meuse.grid
plot(meuse.grid, col="yellow")

data(meuse)
coordinates(meuse) <- c("x", "y")
points(meuse$x, meuse$y)
Benutzeravatar
EDi
Beiträge: 1599
Registriert: Sa Okt 08, 2016 3:39 pm

Re: Grid zuschneiden

Beitrag von EDi »

Dann eine convex hull nehmen (das Polygon, dass durch die äußeren Punkte aufgespannt wird) und dann zuschneiden.

Code: Alles auswählen

library(raster)
library(tidyverse)

# an empty raster of global extent
r <- raster()
# make up some raster values
values(r) <- runif(length(r))

# make up some random coordinates around north america
coords <- cbind( 
  lon = runif(100, min = -120, max = -60),
  lat = runif(100, min = 30, max = 50))

# let's have a look
plot(r)
points(coords, add = TRUE)

# get the convex hull
hull_points <- coords[chull(coords),]

# convert to polygon
hull_polygon <- hull_points %>%
  Polygon() %>%
  list() %>%
  Polygons(1) %>%
  list() %>%
  SpatialPolygons()

# mask the raster
rr <- mask(r, hull_polygon)

# let's have another look
plot(rr)
Quelle:
https://stackoverflow.com/questions/452 ... onvex-hull
Bitte immer ein reproduzierbares Minimalbeispiel angeben. Meinungen gehören mir und geben nicht die meines Brötchengebers wieder.

Dieser Beitrag ist lizensiert unter einer CC BY 4.0 Lizenz
Bild.
Everlast

Re: Grid zuschneiden

Beitrag von Everlast »

Genial, genau das was ich gesucht hab :D danke EDi, auch für den Link!
Antworten