在 R 中使用不规则间隔的网格 netcdf 数据

问题描述

我很快就会处理模拟的空气污染数据。我得到了一些示例数据，看起来网格非常不规则（'brick' 会抛出关于不规则网格的错误），这使得将其转换为栅格格式变得更加困难。

    # load packages ####
    library(ncdf4)
    library(raster)
    library(rgdal)
    library(fields)
    
    # netcdf to data matrix ####
    nc_data <- nc_open('filename.nc')

    cnc_PM2_5       <- ncvar_get(nc_data,attributes(nc_data$var)$names[1])
    longitude       <- ncvar_get(nc_data,"longitude")
    latitude        <- ncvar_get(nc_data,"latitude")
    #time            <- ncvar_get(nc_data,"time")

    # plot matrix,latitude is upside down ####
    image.plot(longitude,rev(latitude),PM2_5,main="PM2_5",ylab="latitude")

    # get quick shapefile for Finland (http://www.diva-gis.org/datadown)####
    finland <- readOGR("FIN_adm4.shp")
    # add map projection
    proj4string(finland) <- CRS("+proj=longlat +datum=wgs84 +init=epsg:4326")
    helsinki <- finland[finland$NAME_4=="Helsinki",]
    # add Helsinki outline to image.plot
    plot(helsinki,add=TRUE)

如何将这些不规则数据转化为空间对象？（SpatialPointsDataFrame 或 raster）

可重现的例子

    PM25 <- matrix(data=c(4.540674,3.437939,3.373072,2.335204,2.140603,2.529804,4.475807,4.346074,6.42181,3.113605,2.594671,2.854138,3.632539,6.097476,3.308205,5.643409,4.151473,2.400071,3.827139,4.994741,3.243339,2.075737,4.281207,3.697406,5.448809,3.567672,2.983871,5.578542,2.20547,6.22721,3.502806,5.773142,5.254208,3.762273,2.724404,3.892006,5.838009,6.746144,5.059608,2.270337,4.086606,4.605541,4.41094,5.319075,6.162343,3.048738,4.02174),ncol=10,nrow=10)
    
    # latitude
    latitude <- c(60.1191177368164,60.1192321777344,60.1193504333496,60.1194686889648,60.1195831298828,60.119701385498,60.1198196411133,60.1199340820312,60.1200523376465,60.1201667785645)
    # longitude
    longitude <- c(24.5949993133545,24.595235824585,24.5954704284668,24.5957050323486,24.5959415435791,24.5961761474609,24.5964126586914,24.5966472625732,24.5968818664551,24.5971183776855)
    
    image.plot(longitude,latitude,PM25,main="PM2_5")

解决方法

感谢@Spacedman 提供了 this 答案：

将存储为矩阵的数据转换为坐标点的向量：

首先将坐标展开成一组完整的坐标点：

latlong = expand.grid(long=longitude,lat=latitude)

然后通过使用 c(.) 展平矩阵来添加包含数据的列：

latlong = cbind(latlong,PM25 = c(PM25))

给予：

head(latlong)
      long      lat PM25
1 24.59500 60.11912    1
2 24.59524 60.11912    2
3 24.59547 60.11912    3
4 24.59571 60.11912    4
5 24.59594 60.11912    5
6 24.59618 60.11912    6

然后用sp设置坐标：

coordinates(latlong)=~long+lat

当绘制时，它应该看起来与 image.plot 一样，除了点，可能还有不同的颜色：

spplot(latlong,"PM25")

geospatial netcdf r r