Como produzir grade espacial a partir de varredura?

9

Preciso de uma grade espacial como uma grade principal para diversos mapas temáticos. Como produzo uma Grade espacial a partir de uma varredura descartando todos os pixels de NA?

Janbob Squarebrains
fonte
6
Jogue-nos alguns pedaços aqui. Um pouco de código para criar uma varredura e o que você espera dela?
Spacedman

Respostas:

14

Você pode obter todas as coordenadas de células não NA em uma varredura com:

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

essas são as células que não são NA. Você provavelmente pode alimentá-los com SpatialPixels

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 

Embora pessoalmente qualquer coisa em uma grade eu manteria como uma varredura.

Ainda não tenho muita certeza do que você deseja - os SpatialGridobjetos definem grades retangulares completas, portanto uma sem os pixels de NA não faz sentido.

Spacedman
fonte
obrigado companheiro, eu também não tenho certeza do que querer. Tento desenvolver um fluxo de trabalho para produzir mapas raster idênticos (em termos de resolução, coorinates, crs, ...) a partir de vários dados de pontos espaciais. rasterize não é uma opção devido à dispersão de pontos espaciais. Precisa usar o iwd ou similar. Você está certo, a grade é retangular - o que eu preciso como modelo mestre provavelmente é um SpatialPointsDataFrame. que pode ser quadriculado. S.
Janbob Squarebrains
Você provavelmente desejará uma varredura mestre (ha!) E criará todas as grades subsequentes simultaneamente a ela.
Radar
2
Acordado. Se as coisas estão em grade, a varredura é quase sempre a resposta. Ou isso, ou percorrendo cerca de 20 páginas da Análise de dados espaciais aplicada no R para descobrir como fazer o SpatialPixels. De fato, eu dei essa escolha a um estudante de doutorado hoje. Ele escolheu ... sabiamente!
Spacedman
20

Para transformar um RasterLayer em um objeto Spatial * (Grid ou Pixels), você pode usar a função de coerção "como"

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
fonte
7

Seus dois requisitos parecem ser sobre coisas diferentes:

1) Algum tipo de modelo de grade raster confiável.

2) Uma grade esparsa que não armazena explicitamente as células ausentes.

O sp :: GridTopology fornece o primeiro, é apenas uma descrição da grade com base na coordenada da célula inferior esquerda (cellcentre.offset), no espaçamento da célula (tamanho da célula) e nas dimensões da grade (cells.dim).

A classe sp :: SpatialPixelsDataFrame permite armazenar grades esparsas, mas, por si só, armazena muito mais que o "modelo" - também armazena todas as coordenadas explicitamente. Isso ocorre porque ele realiza dois trabalhos, um que permite preservar as coordenadas originais que vêm da grade e são possivelmente levemente irregulares; dois permite armazenar apenas as células que possuem valores válidos. (Indiscutivelmente * esses dois objetivos deveriam ter sido separados, mas isso é outra história).

Eu não acho que o pacote raster tenha um análogo explícito ao GridTopology, mas você pode se apossar dos componentes para "rodar sozinho":

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

Vinculando tudo isso, podemos passar de raster para sp:

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

(Observe como as dimensões devem ser revertidas). Outra maneira mais simples é coagir no SpatialGrid e usar o getGridTopology de sp, embora isso seja mais caro, pois as coordenadas acabam sendo geradas ao longo do caminho:

getGridTopology(as(r1, "SpatialGrid"))

Essas três partes da "topologia" da varredura não são todas necessárias, pois algumas são redundantes, mas não existe uma maneira de capturá-las como um único objeto - exceto que a varredura criada acima é "vazia" e, portanto, pode fazer o trabalho que GridTopology faz para sp. Não tenho certeza dos detalhes de quão "vazio" é, mas certamente não armazena explicitamente o slot "data" e é menor do que seria com valores nele. O pacote raster geralmente faz o possível para manter o uso da memória no mínimo, e com isso você pode não precisar se preocupar em ser "escasso".

Isso pode ajudar a explicar um pouco mais, eu sei que estou sobrepondo a resposta de Spacedman, mas ainda não está claro exatamente o que você quer dizer com a pergunta.

  • (Pode-se argumentar que, como você pode armazenar células esparsas apenas armazenando o índice em vez de coordenadas explícitas, as coordenadas de células "ligeiramente irregulares" originalmente poderiam ser armazenadas como atributos se você realmente as quisesse. Nem sp nem raster lidam com rasters irregulares, nem mesmo o caso retilíneo simples - isso está de acordo com a maioria das ferramentas GIS, mas é lamentável, pois é bastante comum em formatos como o NetCDF e é tratado pela boa e velha função graphics :: image de R (embora não se a opção rasterImage recente for usada) .)
mdsumner
fonte
Obrigado pelas declarações muito instrutivas que preciso digerir. Espero que me permita resolver meu problema ou postar uma pergunta mais explícita!
Janbob Squarebrains
Também @Spacedman: eu defini com mais precisão meu histórico de perguntas aqui no fórum.
Janbob Squarebrains