R Studio中统计遥感图像像素数及分割农田影像的技术求助

一、统计遥感图像像素数量 #

无需逐个波段转数据框处理，raster包提供了更高效的统计方法：

library(raster)

# 读取多波段影像（推荐用brick适配多波段数据）
grey <- brick("SentFurrowField2021WB.tif")

# 1. 获取影像总像素数
total_pixels <- ncell(grey)
cat("总像素数：", total_pixels, "\n")

# 2. 统计单波段有效像素（非NA值，以波段1为例）
band1_valid <- sum(!is.na(values(grey[[1]])))
cat("波段1有效像素数：", band1_valid, "\n")

# 3. 统计所有波段共同有效的像素（所有波段均无NA的像素）
all_bands_valid <- sum(complete.cases(values(grey)))
cat("全波段有效像素数：", all_bands_valid, "\n")

二、分割影像为上/中/下三部分并计算NDVI #

1. 影像分割 #

通过影像y轴范围等距划分，用crop函数裁剪出三个区域：

# 获取影像y轴的最小/最大值
y_min <- ymin(grey)
y_max <- ymax(grey)

# 计算三等分阈值
split_upper <- y_min + 2*(y_max - y_min)/3
split_middle <- y_min + (y_max - y_min)/3

# 裁剪三个区域
upper <- crop(grey, extent(xmin(grey), xmax(grey), split_upper, y_max))
middle <- crop(grey, extent(xmin(grey), xmax(grey), split_middle, split_upper))
lower <- crop(grey, extent(xmin(grey), xmax(grey), y_min, split_middle))

2. 计算NDVI #

NDVI公式为(近红外波段 - 红波段)/(近红外波段 + 红波段)，请根据你的影像波段对应关系调整索引（示例假设波段4为近红外，波段3为红波段）：

# 定义NDVI计算函数
calc_ndvi <- function(raster_obj) {
  nir <- raster_obj[[4]]
  red <- raster_obj[[3]]
  ndvi <- (nir - red) / (nir + red)
  return(ndvi)
}

# 对三个区域分别计算NDVI
upper_ndvi <- calc_ndvi(upper)
middle_ndvi <- calc_ndvi(middle)
lower_ndvi <- calc_ndvi(lower)

# 可视化结果
plot(upper_ndvi, main = "上部区域NDVI")
plot(middle_ndvi, main = "中部区域NDVI")
plot(lower_ndvi, main = "下部区域NDVI")

补充说明 #

• 若影像波段对应关系不同，需修改函数中nir和red的波段索引；
• 若需按田块实际边界分割，可使用sf包绘制多边形，再用mask函数提取对应区域。

内容的提问来源于stack exchange，提问作者Grecia Magdaleno