使用原始网格大小从 R 中的 netCDF 文件制作地图

Question

我想使用 ggplot 在 R 中映射 netCDF 文件中的一些数据。

我一直在四处寻找，但没有找到可行的答案。例如，您可以看到这里。

我想使用 netCDF 文件制作地图，该文件保留来自源的预期网格数据（全局数据的~1x1 度网格）。但是，我在使用 ggplot 绘图时遇到了问题，最终需要使用 geom_point() 进行绘图，这使得网格的大小不正确。

我将从上面链接的示例中获取数据。虽然我要稍微调整一下代码，因为我一直在使用“ncdf4”包。

下面是链接示例的改编，其中文件“air.1999.nc”是从同一来源下载的这里。

#Example
library(tidyverse)
library(ncdf4)
library(here)
library(fields)

#Grab the values
temp.nc <- nc_open(here("air.1999.nc"))
temp <- ncvar_get(temp.nc,"air")
ilon <- ncvar_get(temp.nc, 'lon')
ilat <- ncvar_get(temp.nc, 'lat')

#Just pull out a slice for a visual
sample <- temp[,,1,1]

#Check dimensions
dim(sample)

#Take a look
image.plot(sample)

好的，这样通常就可以了，您可以看到预期的网格形式的数据。您可以叠加地图或类似于链接示例中答案的内容。

但是，我想让地图看起来有点不同，并玩转视觉效果。所以我想使用ggplot。

所以我在 YouTube 上找到了一个有趣的解决方案。该代码位于下面，其中包含 YouTube 频道的链接，我在其中获得了下面修改过的代码。

#Function from: https://www.youtube.com/watch?v=OqcYTdSKNYg&t=495s
mapCDF <- function(lat, lon, idat) 
{
  
  #Create a df to plot
  expand.grid(lon, lat) %>%
    dplyr::rename(lon = Var1, lat = Var2) %>%
    mutate(lon = ifelse(lon > 180, -(360 - lon), lon),
           idat = as.vector(idat)) %>% 
    
    #Start plot, feeding in the previous df from expand.grid
    ggplot()+
    geom_point(aes(x = lon, y = lat, color = idat), shape = 'square')+
    scale_color_continuous(name = 'Test Map', na.value = 'transparent')+
  
    #Add the map layer
    borders('world', colour = 'grey15', fill = 'NA')+
    
    #Basic theme settings for background color, text size etc. 
    theme_bw()+
    
    #Map settings, adjusting the lat/long and projection
    coord_map(xlim = c(-170, 170), ylim = c(39, 75), projection = 'orthographic')+ 

    #Remove labels
    ggtitle('')+
    xlab('')+
    ylab('')
  
}

mapCDF(ilat, ilon, sample)

我想你可能会看到这个问题。我已将点的形状从 geom_point() 调整为“正方形”以模仿网格，但它们的大小和间距不正确。他们应该完全填满地图，但他们没有。从默认大小，我得到了同心环效果。如果你将大小调整为，例如，大小 = 5，它会填充更多一点，但这不是我想要的。

可能，geom_point不是正确使用的函数，但我没有运气使用geom_raster或类似的东西，如这个例子。

附加说明：我使用了正交投影，因为我使用的实际数据仅适用于北半球。

非常感谢任何建议。

Answer 1

你是对的，

geom_point

不太适合这个目的。另一方面，

geom_raster

不适用于非笛卡尔投影（因为它将所有点绘制为相同的大小和形状），因此不适用于您的情况。

geom_tile

做我认为你正在寻找的事情。将

geom_point(aes(x = lon, y = lat, color = idat), shape = 'square')+

替换为：

geom_tile(
  aes(x = lon, y = lat, fill = idat, color = idat),
  linewidth = 1
) +

color = idat

美学和

linewidth = 1

论点是在图块之间显示的空白的几个像素上进行绘制。