I would like to map some data from a netCDF file in R using ggplot.

I’ve been looking around and haven’t found a viable answer. You can see, for example, here.

I want to make a map using a netCDF file that retains the anticipated gridded data from the source (~1×1 degree grids for global data). However, I’m running into an issue plotting using ggplot, where I end up needing to use geom_point() to plot, which makes the sizing of the grids incorrect.

I’ll take the data from my linked example above. Though I’m going to adjust the code a bit because I’ve been using the “ncdf4” package.

Below is the adaptation of the linked example, with the file “air.1999.nc” downloaded from the same source here.

#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)

Ok so that generally works and you can see the data in the anticipated grid form. You could superimpose a map or something similar to what the answer had from the linked example.

However, I want to make the maps look a bit different and play around with the visual. So I want to use ggplot.

So I found an interesting solution on YouTube. That code is below with the link to the YouTube channel where I got the code that I have modified below.

#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)

I think you can probably see the issue. I’ve adjusted the shape of the points from geom_point() to “squares” to mimic the grid, but they are sized and spaced incorrectly. They should fill in the map completely, but they don’t. From the default size, I get a concentric ring effect. If you adjust the size to, say, size = 5, it fills in a bit more, but that’s not what I’m looking for.

Probably, geom_point is not the correct function to use, but I’ve had no luck using geom_raster or anything similar, as in this example.

An additional note: I used an orthographic projection because the actual data that I’m using is only for the Northern Hemisphere.

Any advice is greatly appreciated.