Почему растр :: обрезка изменяет значения моего RasterLayer?

Question

У меня возникла любопытная проблема с функцией raster :: crop.Когда я обрезаю большой растр, содержащий уникальные значения идентификатора (равные номеру ячейки), полученный меньший растр содержит дублированные значения, где каждая ячейка должна быть уникальной.

У меня большой растр (r), обрезанный до фрейма SpatialPolygonsData на африканском континенте.Я хочу сгенерировать уникальный идентификатор для каждой ячейки сетки, а затем обрезать больший RasterLayer в меньший многоугольник (одна страна) - таким образом, я смогу впоследствии объединить данные из расчетов меньшего растра в больший растр на основе их уникальной ячейкиидентификаторы.

Однако, когда я обрезаю больший растр, содержащий уникальные идентификаторы, новый растровый слой не содержит уникальных идентификаторов - скорее, существует большое количество дублированных значений.Из того, что я могу сказать, кажется, что некоторые значения ID заменяются соседними значениями, создавая таким образом шаблоны, такие как, например, 301,301,303,303,306,306,306 вместо непрерывного увеличения на 1 (в пределах одной и той же строки курса).Больший r содержит только уникальные значения, поэтому проблема возникает только после обрезки растра в меньший растр.

Проекция одинакова для полигона и растра, и я использую последнюю версию пакета Raster.

Проблема, похоже, возникает только для растров высокого разрешения.Использование вместо этого raster :: mask для замены ячеек вне многоугольника на NA также приводит к дублированию в некоторых случаях.

Проблема для меня полная загадка - я не смог найти никакой возможной причины для этого.Кто-нибудь знает, есть ли проблемы с функцией кадрирования и как она обрабатывает значения?Даже если я найду другой способ сделать это, мне будет интересно, если эта проблема возникает для других растров и, таким образом, ваши данные каким-то образом повреждены в функции обрезки.Я надеюсь, что кто-то может помочь мне выяснить, в чем проблема.

Я создал небольшой пример ниже, который воспроизводит проблему.Если вам нужна дополнительная информация, дайте мне знать.

pack <- c("sp", "raster", "rgdal", "dplyr")
lapply(pack, require, character.only = TRUE); rm(pack)

r <- raster(africa, resolution = 1/60/2) ## Create empty raster layer based on extent of Africa polygons and resolution 30 arc seconds. 
values(r) <- 1:ncell(r) ## Generate unique cell ID (equal to cell number)
poly <- africa[3,] ## Subset one country

r_id <- crop(r, poly) ## Crop r to poly
## This is the function that seems to be responsible for the unexpected result. It should return a smaller raster containing the same values in the same cells as the larger raster. Therefore each grid cell value in r_id should be unique. 

as.data.frame(r_id) %>% ## This just to show that the resulting raster contains duplicate values where none should exist
  group_by(layer) %>% 
  summarise(n = n()) %>% 
  arrange(desc(n))

# A tibble: 137,270 x 2
      layer     n
      <dbl> <int>
 1 24774556     3
 2 24774560     3
 3 24774564     3
 4 24774568     3
 5 24774572     3
 6 24774576     3
 7 24774580     3
 8 24774584     3
 9 24774588     3
10 24774592     3
# ... with 137,260 more rows. 

Answer 1

Вероятно, это связано с неточностью числового значения для больших чисел, когда значения записываются на диск с использованием 4-байтовых чисел с плавающей запятой.

Я подозреваю, что это происходит в вашем случае, потому что r поддерживаетсявверх (временным) файлом.Посмотрите на источник данных после установки номеров ячеек

values(r) <- 1:ncell(r)
r

Я проиллюстрирую это здесь

library(raster)
a <- 24774556:24774565
r <- raster(ncol=2, nrow=5)
values(r) <- a
x <- writeRaster(r, "test1.grd", overwrite=TRUE)
v <- values(x)
v
#[1] 24774556 24774556 24774558 24774560 24774560 24774560 24774562 24774564 24774564 24774564
table(v)
#v
#24774556 24774558 24774560 24774562 24774564 
#       2        1        3        1        3 

Однако, если вы установите тип данных "FLT8S" или "INT4U":

y <- writeRaster(r, "test2.grd", datatype="FLT8S", overwrite=TRUE)
values(y)
#[1] 24774556 24774557 24774558 24774559 24774560 24774561 24774562 24774563 24774564 24774565


z <- writeRaster(r, "test3.grd", datatype="INT4U", overwrite=TRUE)
values(z)
# [1] 24774556 24774557 24774558 24774559 24774560 24774561 24774562 24774563 24774564 24774565

В вашем случае вы можете вместо values(r) <- 1:ncell(r) сделать

 r <- init(r, "cell", datatype="FLT8S", filename="africa_id.grd")

Или пропустить этот бит все вместе.Вы также можете обрезать детали, обрабатывать их, а затем использовать merge, чтобы собрать обрезанные части вместе.