gpt4 book ai didi

r - 如何用不同程度的栅格创建栅格砖?

转载 作者:行者123 更新时间:2023-12-04 03:51:58 25 4
gpt4 key购买 nike

我是的新手R ,所以这个问题很基本,但是我一直在努力,无法找到有效的解决方案。我想从同一区域的一些landat图像创建一个栅格砖。
它们以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

    最佳答案

    以下是一些可以帮助您的事情。由于我没有您的.tif文件,因此仅提供一些提示。您是否检查了栅格的范围?这就是世界的大小,只有这些列,其单元格非常大。因此,在重新采样之前,必须在栅格中添加一个范围。从您的信息中,我做了这样的事情:

    # create an extent that includes all your data
    e<-extent(-46, -43, -2, -0.6)

    # create a raster with that extent, and the number of rows and colums to achive a
    # similar resolution as you had before, you might have to do some math here....
    # as crs, use the same crs as in your rasters before, from the crs slot
    s<-raster(e, nrows=7000, ncols=7800, crs=r1@crs)

    # use this raster to reproject your original raster (since your using the same crs,
    # resample should work fine
    r1<-resample(r1, s, method="ngb")

    节日快乐,


    PS是找到所需范围和分辨率的更好方法:
    # dummy extent from your rasters, instead use lapply(raster list, extent)
    a<-extent(-45.85728, -43.76855, -2.388705, -0.5181549)
    b<-extent(-45.87077, -43.78204, -2.388727, -0.5208711)
    c<-extent(-45.81952 ,-43.7173 , -2.405129 ,-0.5154312)
    extent_list<-list(a, b, c)

    # make a matrix out of it, each column represents a raster, rows the values
    extent_list<-lapply(extent_list, as.matrix)
    matrix_extent<-matrix(unlist(extent_list), ncol=length(extent_list)
    rownames(matrix_extent)<-c("xmin", "ymin", "xmax", "ymax")

    # create an extent with the extrem values of your extent
    best_extent<-extent(min(matrix_extent[1,]), max(matrix_extent[3,]),
    min(matrix_extent[2,]), max(matrix_extent[4,]))

    # the range of your extent in degrees
    ranges<-apply(as.matrix(best_extent), 1, diff)
    # the resolution of your raster (pick one) or add a desired resolution
    reso<-res(r1)
    # deviding the range by your desired resolution gives you the number of rows and columns
    nrow_ncol<-ranges/reso

    # create your raster with the following
    s<-raster(best_extent, nrows=nrow_ncol[2], ncols=nrow_ncol[1], crs=r1@crs)

    关于r - 如何用不同程度的栅格创建栅格砖?,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/20733555/

    25 4 0
    Copyright 2021 - 2024 cfsdn All Rights Reserved 蜀ICP备2022000587号
    广告合作:1813099741@qq.com 6ren.com