r - 绘制选定区域的 netcdf 数据

Question

我发现很难使用我拥有的 netcdf 数据。我浏览了这里的示例(Plotting netcdf file with levels in R)，但我认为我遗漏了一些东西。

我正在尝试绘制 35 度以南区域的混合层深度。数据可以在这里找到(页面底部的最后一个文件):http://www.ifremer.fr/cerweb/deboyer/mld/Surface_Mixed_Layer_Depth.php

该文件有 7 个变量，每个变量包含纬度、经度、时间(12 个月)以及混合层深度的值。

到目前为止我已经:

    MLD <- "mld_DReqDTm02_c1m_reg2.0.nc"
    MLD <- nc_open(MLD)
    print(MLD)

有 7 个变量，我只想要“mld”

lon <- ncvar_get(MLD, varid = "lon")
lat <- ncvar_get(MLD, varid = "lat")

summary(lon)
summary(lat)



MLD$dim$time$units

MLD_1.array <- ncvar_get(MLD, "mld")
dim(MLD_1.array)

length(lon)
length(lat)

ndvi.slice <- MLD_1.array[, , 12] 
dim(ndvi.slice)

mld.vec.long <- as.vector(MLD_1.array)
length(mld.vec.long)

nlon <- dim(lon)
nlat <- dim(lat)
lonlat <- expand.grid(lon, lat)

t <- ncvar_get(MLD, "time")
tunits <- ncatt_get(MLD, "time", "units")
nt <- dim(t)

dname <- "mld" 

tmp.mat <- matrix(mld.vec.long, nrow = nlon * nlat,    ncol = nt)
dim(tmp.mat)

head(na.omit(tmp.mat))

lonlat <- expand.grid(lon, lat)
tmp.df02 <- data.frame(cbind(lonlat, tmp.mat))
 names(tmp.df02) <- c("lon", "lat", "Jan", "Feb", "Mar", "Apr", "May", 
"Jun", "Jul", "Aug", "Sep", "Oct", "Nov", "Dec")

library(reshape)
tmp.df03 <- melt(tmp.df02, id=c("lat","lon"))

tmp.df04 <- subset(tmp.df03, lat >= "-35" & lat <="-80")

 tmp.df04[tmp.df04 ==1.000e+09] <- NA
 tmp.df04
 summary(tmp.df04)


 target = c(-180, 180, -90, -20)
 w <- crop(rgeos::gBuffer(spTransform(countriesLow, CRS(pprj)), width = 0))

我可以在将数据提取到 .csv 后使用 ggplot 绘制此图，删除掩码值，重新导入它然后重新绘制它(需要一段时间......)。

有没有一种方法可以提取一个变量(mld)，然后找到每个月(时间)的最大值并绘制出来？我知道我的代码真的很困惑而且断断续续......

任何帮助将不胜感激!谢谢你!

Answer 1

希望这有帮助!

library(ncdf4)
library(raster)

#Reading netcdf file and extracting lat, lon and variable
MLD <- nc_open("D:/Personal/test/mld_DReqDTm02_c1m_reg2.0.nc")
lon <- ncvar_get(MLD, varid = "lon")
lat <- ncvar_get(MLD, varid = "lat")
mld <- ncvar_get(MLD, "mld")

# using 'raster' package to read the temporal values of variable into a raster
#r1<-flip(raster(t(matrix(mld[,,1], nrow = 180,    ncol = 90))),direction="y")
e<-extent(min(lon),max(lon),min(lat),max(lat))
#extent(r1)<-e
R<-stack()
# creating raster stack of the time series data
for(i in 1:12){
  r1<-flip(raster(t(matrix(mld[,,i], nrow = 180,    ncol = 90))),direction="y")
  extent(r1)<-e
  R<-stack(R,r1)
} # IMPORTANT: This raster stack could be additionally cropped to extract the user's area of interest.
plot(R)

# Extracting max/min values for each time (raster layer) into a dataframe.
df<-data.frame(Months=month.abb)
df$Months <- factor(df$Months, levels = df$Months)
df$months<-c(1:12)
df$MLD_max<-maxValue(R)
df$MLD_min<-minValue(R)
# using 'ggplot2' package
ggplot(df)+geom_point(aes(Months,MLD_min))+geom_line(aes(months,MLD_min))