Лучший способ прорисовать слои растрового кирпича / стека - PullRequest
Новая
       38

Лучший способ прорисовать слои растрового кирпича / стека

0 голосов
/ 01 октября 2018

Я пытаюсь построить все значения в растровом кирпиче для определенных точек.Это делается для создания спектральных графиков для данных дистанционного зондирования для конкретных пикселей.

Я могу сделать это разными способами, но они довольно громоздкие и медленные (см. Пример ниже).Это происходит медленно, потому что преобразование больших растровых файлов в матрицу требует большого объема памяти.

Есть ли лучший способ сделать это с помощью baseR или аккуратного стиха или встроенного способа сделать это в одном из растровых изображений?/ пакеты дистанционного зондирования?

Вот воспроизводимый пример: 

library (raster)
library (rgdal)
library (sp)
library (tidyr)
library (ggplot2)
library (dplyr)


##############################

### dummy data

##############################


coord1 <- c(50, 47, 45)
coord2 <- c(50, 51, 49)
frame <- data.frame(coord1, coord2)
coords <- coordinates(frame)

x = list(1, 2, 3, 4, 5)
y = list(1, 2, 3, 4, 5)

for (i in 1:length(x)) { # create clusters based on coords

set.seed(i+100)
x[[i]] <- rnorm(5000, mean = coords[, 1], sd = 1)
y[[i]] <- rnorm(5000, mean = coords[, 2], sd = 1)
}

obs <- data.frame(x, y)
names(obs) <- c('x1', 'x2', 'x3', 'x4', 'x5', 'y1', 'y2', 'y3', 'y4', 'y5')
coordinates(obs) <- c('x1', 'y1') # this creates a spatial points data frame

# create blank raster of study area

ex <- extent(c(45, 50, 48, 52))
target <- raster(ex, resolution = 0.5)

# create raster brick

r_list <- list()

for (i in 1:ncol(obs)) {

   r_list[[i]] <- rasterize(obs, target, fun = 'count')
}

obs_frequency <- brick(r_list)

И вот одно из возможных, но медленных решений 

############################

### Example Solution

############################

vals_all <- obs_frequency[, ] ### this is very slow ###
vals_all <- data.frame(vals_all)

### sample values

points <- sample_n(vals_all, size = 5)
points <- t(points)
points <- data.frame(points)
points_tidy <- gather(points)
points_tidy$xval <- rep(1:8, times = 5)


### plot

ggplot(points_tidy, aes(x = xval, y = value)) + geom_line(aes(colour = key)) + theme_classic()

1 Ответ

0 голосов
/ 01 октября 2018

Я нашел лучшее решение для этого, используя функцию raster :: extract.Это непосредственно определяет значения и позволяет избежать превращения всего растрового кирпича в матрицу, разрушающую память.

Стоит отметить, что здесь использование кирпича НАМНОГО быстрее, чем использование растрового стека.

############################

### Extract values and plot 

############################

### extract values

points <- c(49, 50, 51) #arbitrary points
pointvals <- raster::extract(obs_frequency, points) ##### USE THE RASTER::EXTRACT FUNCTION

### manipulate data structure

pointvals <- data.frame(t(pointvals))
points_tidy <- gather(pointvals)
points_tidy$xval <- rep(1:8, times = 3)

### plot

ggplot(points_tidy, aes(x = xval, y = value)) + geom_line(aes(colour = key)) + theme_classic()
...