我正在使用卫星数据来观察海洋表面的光照水平。 (来自here的任何数据都有效,因为我只是想了解如何打开和使用数据)
在下面的代码中,我将所有
.nc
文件保存在桌面上的一个文件夹中,然后将其全部导入到 files
下,然后启动循环来打开文件。
有没有办法 1) 裁剪数据
extent(-180, -176, -40, -36)
然后 2) 提取纬度、经度和面值?我猜测有多个深度可用,但我不知道如何检查。
library(raster)
setwd("/Users/AquaModis_PAR")
nc.data<-list()
files<-list.files(pattern = ".nc",full.names=TRUE)
df=NULL
for(i in seq_along(files)) {
par[i] = brick(files[i], stopIfNotEqualSpaced = FALSE, varname = "par")
}
要加载和裁剪数据,请执行以下操作:
library(raster)
files <- list.files(pattern = ".nc", full.names=TRUE)
region <- extent(-180, -176, -40, -36)
data <- lapply(files, function(path) {
# Load and crop.
d <- brick(path, stopIfNotEqualSpaced = FALSE, varname = "chlor_a") |>
crop(region)
# Get raster values.
values <- getValues(d)
# Get coordinates.
coords <- coordinates(d)
cbind(coords, values) |>
data.frame() |>
setNames(c("lon", "lat", "value")) |>
na.omit()
})
查看第一个文件中的数据:
> head(data[[1]])
lon lat value
28 -178.8542 -36.02084 0.06653178
29 -178.8125 -36.02084 0.06850998
30 -178.7708 -36.02084 0.07330192
31 -178.7292 -36.02084 0.07878597
34 -178.6042 -36.02084 0.07993947
35 -178.5625 -36.02084 0.08084643
可能有更优雅的方法将数据提取到数据框中,但这可以完成工作。
最后的
na.omit()
是删除缺少值的记录(我查看的数据中有不少)。
为了完整起见,这些是我测试过的文件:
AQUA_MODIS.20030101.L3m.DAY.CHL.chlor_a.4km.nc
AQUA_MODIS.20030102.L3m.DAY.CHL.chlor_a.4km.nc
AQUA_MODIS.20030103.L3m.DAY.CHL.chlor_a.4km.nc