Я новичок в R, поэтому это очень простой вопрос, но я боролся с ним и не мог найти решение, которое сработало. Я хочу создать растровый кирпич из нескольких изображений Landsat той же области. Они были загружены в формате HDF-EOS
, и я использовал Modis Reprojection Tool, чтобы преобразовать их в .tif
.
Полученные растры имеют одинаковую проекцию, но различаются по размеру, разрешению и происхождению.
Я пробовал несколько подходов, кратко изложенных ниже:
определение экстента подмножества вручную и подмножество всех растров. Затем пытаемся сделать кирпич из подмножества растров.
Повторная выборка растров, чтобы получить в них одинаковое количество столбцов и строк. В идеале это обеспечит выравнивание ячеек растра и их можно поместить в растровый кирпич. Эта опция создавала кирпич, в котором растры не имели значений, они были пустыми.
Мне интересно, какой концепции я должен следовать, чтобы исправить степень. Было бы правильно (и эффективно) создать пустой растр, который я позже заполню значениями импортированного изображения Landsat? Вы видите, где я делаю ошибку? Если это актуально, я работаю над Mac OSX версии 10.9.1 и использую rgdal версии 0.8-14
Любая помощь будет очень признательна!
Спасибо
Я добавляю сюда код, который использовал:
# .tif files have been creating using the Modis Reprojection Tool. Input
# files used for this Tool was LANDSAT HDF-EOS imagery.
library(raster)
library(rgdal)
setwd()=getwd()
# Download the files from dropbox:
dl_from_dropbox <- function(x, key) {
require(RCurl)
bin <- getBinaryURL(paste0("https://dl.dropboxusercontent.com/s/", key, "/", x),
ssl.verifypeer = FALSE)
con <- file(x, open = "wb")
writeBin(bin, con)
close(con)
message(noquote(paste(x, "read into", getwd())))
}
dl_from_dropbox("lndsr.LT52210611985245CUB00-vi.NDVI.tif", "qb1bap9rghwivwy")
dl_from_dropbox("lndsr.LT52210611985309CUB00-vi.NDVI.tif", "sbhcffotirwnnc6")
dl_from_dropbox("lndsr.LT52210611987283CUB00-vi.NDVI.tif", "2zrkoo00ngigfzm")
# Create three rasters
tif1 <- "lndsr.LT52210611985245CUB00-vi.NDVI.tif"
tif2 <- "lndsr.LT52210611985309CUB00-vi.NDVI.tif"
tif3 <- "lndsr.LT52210611987283CUB00-vi.NDVI.tif"
r1 <- raster(tif1, values=TRUE)
r2 <- raster(tif2, band=1, values=TRUE)
r3 <- raster(tif3, band=1, values=TRUE)
### Display their properties
# projection is identical for the three rasters
projection(r1)
# "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
projection(r2)
projection(r3)
# Extents are different
extent(r1)
# class : Extent
# xmin : -45.85728
# xmax : -43.76855
# ymin : -2.388705
# ymax : -0.5181549
extent(r2)
# class : Extent
# xmin : -45.87077
# xmax : -43.78204
# ymin : -2.388727
# ymax : -0.5208711
extent(r3)
# class : Extent
# xmin : -45.81952
# xmax : -43.7173
# ymin : -2.405129
# ymax : -0.5154312
# origin differs for all
origin(r1)
# 5.644590e-05 -8.588605e-05
origin(r2)
# 0.0001122091 -0.0001045107
origin(r3)
# 6.949976e-05 -5.895945e-05
# resolution differs for r2
res(r1)
# 0.0002696872 0.0002696872
res(r2)
# 0.0002696875 0.0002696875
res(r3)
# 0.0002696872 0.0002696872
## Try different approaches to create a raster brick
# a- define a subset extent, and subset all the rasters
plot(r1, main="layer1 NDVI")
de <- drawExtent(show=TRUE, col="red")
de
# class : Extent
# xmin : -45.36159
# xmax : -45.30108
# ymin : -2.002435
# ymax : -1.949501
e <- extent(-45.36159,-45.30108,-2.002435,-1.949501)
# Crop each raster with this extent
r1c <- crop(r1,e)
r2c <- crop(r2,e)
r3c <- crop(r3,e)
# Make raster brick
rb_a <- brick(r1c,r2c,r3c)
# Error in compareRaster(x) : different extent
# b- Resample each raster
s <- raster(nrow=6926, ncol=7735) # smallest nrow and ncol among r1,r2 and r3
r1_res <- resample(r1,s, method="ngb")
r2_res <- resample(r2,s, method="ngb")
r3_res <- resample(r3,s, method="ngb")
# Resampling gives for the three rasters the following message:
# Warning message:
# In .local(x, y, ...) :
# you are resampling y a raster with a much larger cell size,
# perhaps you should use "aggregate" first
# Make raster brick
rb_c <- brick(r1, r2, r3)
# Error in compareRaster(x) : different extent
library(RCurl)
? - person thijs van den bergh   schedule 23.12.2013