我的目标是在世界地图上绘制硝酸盐 (no3) 数据,使用这些数据的正确经度和纬度。
有两个 netcdf 文件:
1. 与数据
2. 与网格信息
数据摘要信息:
no3 是长度为 x*y*sigma 的数组
no3_df 是 'x*y 观测值。 3 个变量'
x = 整数 [180]
y = 整数 [193]
西格玛 = 数组[53]
我想查看 sigma(“深度”)20。因此,我执行了以下操作:
# Load the needed libraries to handle netcdf files
library(ncdf)
library(akima)
# Open data and grid files
file1 <- open.ncdf(file.choose())
grid <- open.ncdf(file.choose())
# Read relevant variables/parameters from data file1
x <- get.var.ncdf(file1,varid="x")
y <- get.var.ncdf(file1,varid="y")
sigma <- get.var.ncdf(file1,varid="sigma")
no3 <- get.var.ncdf(file1,varid="no3")
sigma_plot <- no3[,,sigma=20]
# Read relevant variables/parameters from grid file
plon <- get.var.ncdf(grid,varid="plon")
plat <- get.var.ncdf(grid,varid="plat")
# Each cell of sigma_plot corresponds to one cell of plon and plat.
A <- array(c(plon,plat,sigma_plot),dim=c(180,193,3))
# Now B is an array containing for each row: (longitude, latitude, value).
B <- apply(A, 3, cbind)
# But it is not a regular grid, so interpolate to a regular grid. akima library
C <- interp(B[,1],B[,2],B[,3],
xo=seq(-180,180,1),yo=seq(-90,90,by=1), # tweak here the resolution
duplicate='mean') # extra y values are duplicates
#########
# PLOTTING
#########
# This one works, but doesn't have a correct longitude and latitude:
filled.contour(x,y,sigma_plot, col=rich.colors(18))
# Try to plot with lon and lat
filled.contour(C, col=rich.colors(30))
由于 fill.contour 图没有正确的经度和纬度,我想使用 ggplot。但是,我不知道该怎么做......
# And the plotting with ggplot
ggplot(aes(x=plon_datafrm,y=plat_datafrm),data=no3_df) +
geom_raster() +
coord_equal() +
scale_fill_gradient()
这似乎不起作用。我对 ggplot 很感兴趣,所以这可能就是原因,我真的很感激任何帮助。