如何将我的时空 NetCDF 数据转换为空间数据？

Question

我是 R 语言的初学者。我完全被这个问题所困扰。您可以从下面的链接下载 netCDF 文件进行查看。

https://drive.google.com/file/d/0ByY3OAw62EShbkF6VWNFUkRYMmM/view?usp=sharing

^这是我的 NetCDF 大气数据文件，有 8 个变量和 8 个维度。在这里，我感兴趣的变量是：

TIMSID 是站点编号（包括城市站点、农村站点等）
URBAN :: 城市站点数 [urban 是 3 行 250 列矩阵。第1行是城市站点数，第2行是纬度，第3行是经度。]
TIME :: 数据收集时间为 2012 年 3 月 1 日至 2012 年 5 月 [“时间”编码为 YYYYMMDDHH]
PM10 :: 每个站点每个站点测量的每小时颗粒物浓度

我只需要处理这个大型数据集中的这 4 个变量。

我必须仅在“2012 年 3 月 1 日”的城市站点中分离 PM10 值的数据。 （实际上我需要在 TIMSID 变量中找到哪些站点是城市站点，并且仅匹配 2012 年 3 月 1 日城市站点的相应 PM10 值。）

例如，在 TIMSID 中，存在不同类型的站点，例如城市、农村等，名称为 111121,111122,111123,111124，但城市站点编号为 111121,111123..等，因此我必须仅考虑 TIMSID 数据中的城市站点和想匹配对应的pm10值，时间，纬度，经度。然后终于想做一个新的数据集。

最终的表格/数据集应该是 ~ column1-time（仅 2012 年 3 月 1 日），column2-城市站点编号，第 (3,4) 列-相应城市站点的纬度和经度，第 5 列-每个小时的 pm10 值城市站点

我已使用以下命令从 NetCDF 文件中读取数据。但是我不明白我应该进一步做什么......

install.packages("ncdf",dependencies=TRUE)
library(ncdf)

nc<-open.ncdf("2012_03_05_PM10_surface.nc")
print(nc)

tmsid<-get.var.ncdf(nc,"TMSID")
timsid

urban<-get.var.ncdf(nc,"urban")
urban
time<-get.var.ncdf(nc,"TIME")

pm10<-get.var.ncdf(nc,"PM10")

由于我是 R 初学者，所以我只知道基本命令。我不知道，我应该学习哪个特定的包来解决这个问题。请帮帮我？提前感谢您的宝贵时间。如果您需要任何进一步的信息，请随时询问我。

Answer 1

library(ncdf)
nc <- open.ncdf("2012_03_05_PM10_surface.nc")
tmsid <- get.var.ncdf(nc,"TMSID")
urban <- get.var.ncdf(nc,"urban")
time <- get.var.ncdf(nc,"TIME")
pm10 <- get.var.ncdf(nc,"PM10")

我们先来看看nc：

[1] "file ~/Downloads/2012_03_05_PM10_surface.nc has 8 dimensions:"
[1] "data_num   Size: 683016"
[1] "ncl1   Size: 683016"
[1] "obsnum_urban   Size: 250"
[1] "ID_LAT_LON   Size: 3"
[1] "obsnum_road   Size: 33"
[1] "obsnum_background   Size: 5"
[1] "obsnum_rural   Size: 16"
[1] "ncl7   Size: 683016"
[1] "------------------------"
[1] "file ~/Downloads/2012_03_05_PM10_surface.nc has 8 variables:"
[1] "int TMSID[data_num]  Longname:TMSID Missval:NA"
[1] "int TIME[ncl1]  Longname:TIME Missval:NA"
[1] "float PM10[data_num]  Longname:PM10 Missval:1e+30"
[1] "float urban[ID_LAT_LON,obsnum_urban]  Longname:urban Missval:1e+30"
[1] "float road[ID_LAT_LON,obsnum_road]  Longname:road Missval:1e+30"
[1] "float background[ID_LAT_LON,obsnum_background]  Longname:background Missval:1e+30"
[1] "float rural[ID_LAT_LON,obsnum_rural]  Longname:rural Missval:1e+30"
[1] "int TMS_JULIAN[ncl7]  Longname:TMS_JULIAN Missval:NA"

它告诉我们的是urban 的行是ID、纬度和经度。然后我们有tmsid 给出与time 的向量相同大小的ID 向量：每个data_num 一个，即。 e. PM10 中每个数据点的一对 ID 时间，这意味着我们将能够通过 ID（由 urban 的第一行给出）和时间戳（从 2012030101 到 2012030124）对 pm10 进行子集化。

# First we need to make a dataframe out of urban, for convenience.
urban <- as.data.frame(t(urban))
colnames(urban) <- c("ID", "LAT", "LON")
# Then we do the subsetting using a lapply, so we can batch-subset:
res <- lapply(urban$ID, 
              function(x)data.frame(ID=x,
                                    pm=pm10[tmsid%in%x & time%in%2012030101:2012030124], 
                                    time=2012030101:2012030124))
# Which gives us a list of sub-dataframes that we want to compress back into a single dataframe:
res <- do.call(rbind,res)
# Finally we merge that with the original urban dataframe
# so that each entry has its own LAT and LON:
res <- merge(res, urban, by="ID")
res
#         ID   pm       time      LAT      LON
#1    111121   42 2012030101 37.56464 126.9760
#2    111121   36 2012030102 37.56464 126.9760
#3    111121   46 2012030103 37.56464 126.9760
#4    111121   40 2012030104 37.56464 126.9760
#5    111121   36 2012030105 37.56464 126.9760
#...
#5995 831154   81 2012030119 37.52662 126.8064
#5996 831154   72 2012030120 37.52662 126.8064
#5997 831154   81 2012030121 37.52662 126.8064
#5998 831154   70 2012030122 37.52662 126.8064
#5999 831154   74 2012030123 37.52662 126.8064
#6000 831154   74 2012030124 37.52662 126.8064

250 个城市站点 X 24 小时 = 6000 个数据点，这确实是我们得到的。