9

我想在不手动定义范围的情况下裁剪一些栅格的无数据部分( 1中无数据为黑色的图像示例)。

任何的想法?

没有数据的图像

4

5 回答 5

21

您可以使用trim删除仅具有NA值的外部行和列:

library(raster)
r <- raster(ncols=18,nrows=18)
r[39:49] <- 1
r[205] <- 6
s <- trim(r) 

要更改其他值或更改其他值,NA您可以使用reclassify. 例如,要更改NA为 0:

x <- reclassify(r, cbind(NA, 0))
于 2015-12-05T03:48:43.730 回答
4

[为对象定义了子集和[<-替换方法,raster因此您可以简单地r[ r[] == 1 ] <- NA摆脱1nodata 值所在的值(NAvalue(r)如果您不确定,请使用找出 R 认为您的 nodata 值应该是什么)。

请注意,您必须r[][subsetting 命令中使用才能访问这些值。这是一个工作示例...

例子

#  Make a raster from system file
logo1 <- raster(system.file("external/rlogo.grd", package="raster"))

#  Copy to see difference
logo2 <- logo1

#  Set all values in logo2 that are > 230 to be NA
logo2[ logo2[] > 230 ] <- NA

#  Observe difference
par( mfrow = c( 1,2 ) )
plot(logo1)
plot(logo2)

在此处输入图像描述

于 2013-10-30T10:50:15.073 回答
4

我有 2 个略有不同的解决方案。第一个需要手动识别范围,但使用预定义的功能。第二个更自动,但更手工制作。

创建前 2 行为 NA 的可重现栅格

library(raster)
# Create a reproducible example
r1 <- raster(ncol=10, nrow=10)
# The first 2 rows are filled with NAs (no value)
r1[] <- c(rep(NA,20),21:100)

解决方案#1

使用 drawExtent() 从绘制的图形中手动获取范围

plot(r1)
r1CropExtent <- drawExtent()

使用从图中选择的范围裁剪栅格

r2 <- crop(r1, r1CropExtent)

作比较

layout(matrix(1:2, nrow=1))
plot(r1)
plot(r2)

解决方案#2

它标识栅格中仅具有 NA 值的行和列,并删除位于栅格边缘的行和列。然后它使用 计算范围extent()

将栅格转换为一个矩阵,用于标识值是否为 NA。

r1NaM <- is.na(as.matrix(r1))

查找未完全被 NA 填充的列和行

colNotNA <- which(colSums(r1NaM) != nrow(r1))
rowNotNA <- which(rowSums(r1NaM) != ncol(r1))

通过使用未完全由 NA 填充的第一个和最后一个列和行来查找新栅格的范围。用于crop()裁剪新栅格。

r3Extent <- extent(r1, rowNotNA[1], rowNotNA[length(rowNotNA)],
   colNotNA[1], colNotNA[length(colNotNA)])
r3 <- crop(r1, r3Extent)

绘制栅格以进行比较。

layout(matrix(1:2, nrow=1))
plot(r1)
plot(r3)
于 2013-10-31T12:15:38.300 回答
1

我已经根据 Marie 的回答编写了一个小函数来快速绘制裁剪的栅格。但是,如果栅格非常大,则可能存在内存问题,因为计算机可能没有足够的 RAM 来将栅格加载为矩阵。

因此,我编写了一个内存安全函数,如果计算机有足够的 RAM(因为它是最快的方法),它将使用 Marie 的方法,或者如果计算机没有足够的 RAM(它速度较慢但内存安全),则使用基于光栅函数的方法)。

这是功能:

plotCroppedRaster <- function(x, na.value = NA)
{
  if(!is.na(na.value))
  {
    x[x == na.value] <- NA
  }
  if(canProcessInMemory(x, n = 2))
  {
    x.matrix <- is.na(as.matrix(x))
    colNotNA <- which(colSums(x.matrix) != nrow(x))
    rowNotNA <- which(rowSums(x.matrix) != ncol(x))

    croppedExtent <- extent(x, 
                            r1 = rowNotNA[1], 
                            r2 = rowNotNA[length(rowNotNA)],
                            c1 = colNotNA[1], 
                            c2 = colNotNA[length(colNotNA)])

    plot(crop(x, croppedExtent))
  } else
  {
    xNA <- is.na(x)
    colNotNA <- which(colSums(xNA) != nrow(x))
    rowNotNA <- which(rowSums(xNA) != ncol(x))

    croppedExtent <- extent(x, 
                            r1 = rowNotNA[1], 
                            r2 = rowNotNA[length(rowNotNA)],
                            c1 = colNotNA[1], 
                            c2 = colNotNA[length(colNotNA)])

    plot(crop(x, croppedExtent))
  }
}

例子 :

library(raster)
r1 <- raster(ncol=10, nrow=10)
r1[] <- c(rep(NA,20),21:100)

# Uncropped
plot(r1)
# Cropped
plotCroppedRaster(r1)

# If the no-data value is different, for example 0
r2 <- raster(ncol=10, nrow=10)
r2[] <- c(rep(0,20),21:100)

# Uncropped
plot(r2)
# Cropped
plotCroppedRaster(r2, na.value = 0)
于 2015-12-02T12:06:46.317 回答
0

如果您使用该rasterVis软件包(2021 年 6 月 25 日之后的任何版本),它将自动裁剪NAterra's的值SpatRaster

rasterVis从 GitHub安装开发版本

if (!require("librarian")) install.packages("librarian")
librarian::shelf(raster, terra, oscarperpinan/rastervis)
# Create a reproducible example
r1 <- raster(ncol = 10, nrow = 10)

# The first 2 rows are filled with NAs (no value)
r1[] <- c(rep(NA, 20), 21:100)

levelplot()为了r1

rasterVis::levelplot(r1,
                     margin = list(axis = TRUE))

转换为terra's SpatRaster然后再次使用levelplot()

r2 <- rast(r1)
rasterVis::levelplot(r2,
                     margin = list(axis = TRUE))

reprex 包于 2021-06-26 创建 (v2.0.0 )

于 2021-06-27T05:53:00.143 回答