如何从 netcdf 文件可视化地图？答案

【问题标题】：How to visualize a map from a netcdf file?如何从 netcdf 文件可视化地图？
【发布时间】：2012-08-30 09:22:00
【问题描述】：

我有一个 netcdf 文件，我想将土壤深度图可视化

   [1] "file C:\\Users\\SoilDepth-gswp.nc has 3 dimensions:"
     [1] "x   Size: 360"
     [1] "y   Size: 150"
     [1] "land   Size: 15238"
     [1] "------------------------"
     [1] "file C:\\SoilDepth-gswp.nc has 3 variables:"
     [1] "float nav_lon[x,y]  Longname:Longitude Missval:1e+30"
     [1] "float nav_lat[x,y]  Longname:Latitude Missval:1e+30"
     [1] "float SoilDepth[land]  Longname:Soil depth Missval:1.00000002004088e+20"

似乎我必须将纬度与经度以及土地点连接起来才能获得土壤深度的地图。我真的很困惑。谁能帮我处理这种数据。

【问题讨论】：

您的网格大小为 (360 * 150 = 54e3)，而您的 land 变量的大小为 15238，这不是您的网格大小的倍数。你对此有什么解释吗？
你从哪里得到这些数据？还是你自己创造的？差异可能是由于存在 NA 值，即海洋中没有土壤深度。
land 不是变量而是另一个维度，因此它的大小不必是 x * y
@plannapus 你是对的。我不完全确定文件中的内容，但我认为这是某个研究人员制作的特定格式的文件，您需要他/她的一些代码才能正确解释文件的内容。我会尝试联系文件的创建者。
投票关闭是因为过于本地化，这是这个特定文件特有的问题。

标签： r dictionary ggplot2 netcdf

【解决方案1】：

我更喜欢使用ggplot2 包进行可视化。使用@plannapus 的优秀解决方案：

require(reshape)
require(ggplot2); theme_set(theme_bw())
land_df = melt(land)
ggplot(aes(x = X1, y = X2, fill = value), data = land_df) + 
  geom_raster() + coord_equal() + 
  scale_fill_continuous(na.value = "transparent")

如果要更改轴的标题，请不要更改aes 中的变量名称。这些值引用数据中的列，更改它们会导致您得到错误，land_df 中没有名为X 的轴。如果要更改放置在轴上的名称：

ggplot(aes(x = X1, y = X2, fill = value), data = land_df) + 
  geom_raster() + coord_equal() + 
  scale_fill_continuous(na.value = "transparent") + 
  scale_x_continuous("X")

【讨论】：

如果它们打扰到您，您应该可以通过以下方式删除它们：+theme(axis.title=element_blank())

【解决方案2】：

library(ncdf)
# I'm assuming this is the netcdf file you are working with:
download.file("http://dods.ipsl.jussieu.fr/gswp/Fixed/SoilDepth.nc", destfile="SoilDepth.nc")
soil <- open.ncdf("SoilDepth.nc")
#The way to extract a variable is as following:
soil$var[[3]] -> var3 # here, as shown in your question, SoilDepth is the 3rd variable
get.var.ncdf(soil, var3) -> SoilDepth

dim(SoilDepth)
[1] 15238

正如您在 netcdf 文件的摘要中所说，变量 SoilDepth 仅取决于维度 land，而不取决于 x 和 y，所以我不确定当它出现时它会给您带来什么影响绘制这个数据集。

编辑

原来有一个键链接x、y和land：

download.file("http://dods.ipsl.jussieu.fr/gswp/Fixed/landmask_gswp.nc", destfile="landmask.nc")
landmask <- open.ncdf("landmask.nc")
landmask$var[[3]] -> varland
get.var.ncdf(landmask, varland) -> land
sum(land==1)
[1] 15238

所以为了绘图：

# The values where stored in an expected order, hence the transpose
land = t(land)
land[land==1] <- SoilDepth
land[land==0] <- NA
land = t(land)
image(land)

【讨论】：

仅在 get.var.ncdf(nc, "SoilDepth") 的结果上调用 plot 也不会产生好的时间序列。奇怪的文件...
维度 land 和维度 x 和 y 之间可能存在隐含关系，但我似乎找不到它......
我也很难找到关系。 OP 需要询问数据的创建者...
好吧，也许这个文件在这里：dods.ipsl.jussieu.fr/gswp/Fixed/landmask_gswp.nc 是关键。
如果您有新问题，请发布新问题。

【解决方案3】：

你想用 R 可视化它吗？

如果使用其他软件进行可视化没有问题，您可以使用 ncBrowse，可用here，或 Panoply，一种更复杂的软件，由 NASA 提供，您可以下载here。

如果你想使用 R，你可以使用 ncdf 包。借助get.var.ncdffunction，您将能够提取数据。您可以通过sp 包和spplotfunction 或使用rglpackage（或scatterplot）来绘制它。

【讨论】：

我在@plannapus 的回答中添加了一些代码来解决这个问题

【解决方案4】：

为了快速查看文件，您可以使用 ncview。这些地图不是特别漂亮，但对于了解给定文件的外观非常有用。这也很容易在远程服务器上运行。

请看这里：http://meteora.ucsd.edu/~pierce/ncview_home_page.html

【讨论】：