CODEWITHSUNDEEP
rconvertersnetcdf
Douie
Douie

Reputation: 211

Convert NetCDF to SpatialGridDataFrame with R

I try to convert a netCDF file to a SpatialGridDataFrame class in order to plot a map. In my netCDF file I have the coordinates (lon/lat), the time (1431 values) and the value of relative humidity (777232 values which corresponds to 532 points (lat/lon grid) and 1431 moments (which represents 4 measures in a day during one year)). I want to plot a map for one moment ie the values of relative humidity on a grid of 532 points.

I don't really know the difference between SpatialGridDataFrame or SpatialPixelsDaraFrame, I tried to obtain the results with both functions. Here is my code:

library(ncdf)
library(maptools)

setwd("~/Documents/Data")


#Lire les données du ncdf, dans le cas d'une seule variable et 3 dimensions - Read data
nc <- open.ncdf("X195.221.112.194.318.6.28.31.nc")
lat <- get.var.ncdf(nc,"lat")
long <- get.var.ncdf(nc,"lon")
time <- get.var.ncdf(nc,"time")
data <- get.var.ncdf(nc)


#Redimensionner les données dans le bon format - To have good dimensions
lonlat <- SpatialPoints(expand.grid(long,lat), proj4string=CRS(as.character(NA)))
data <- as.data.frame(matrix(data,ncol=1461))

#Pas de temps pour lequel on veut tracer la carte - select time for the map
i=1461

#Changer la classe de l'objet de NetCDF vers SP - convert to SP
Grid <- SpatialPixelsDataFrame(lonlat,data[[i]], tolerance = sqrt(.Machine$double.eps))

And this is the error I got:

invalid class “SpatialPixelsDataFrame” object: invalid object for slot "data" in class "SpatialPixelsDataFrame": got class "numeric", should be or extend class "data.frame"
Warning:
In points2grid(points, tolerance, round) :
  grid has empty column/rows in dimension 1

An idea to help?

After the answer:

library(ncdf)
library(maps)
library(fields)
library(scales)

setwd("~/Documents/Data")


nc <- open.ncdf("X195.221.112.194.322.0.48.19.nc")
lat <- get.var.ncdf(nc,"lat")
long <- get.var.ncdf(nc,"lon")
time <- get.var.ncdf(nc,"time")
var <- get.var.ncdf(nc, "hgt")
close.ncdf(nc)

i=100
nom_var="Hgt"
var.slice <- var[,,i]

nlevel=64   #number of colors
image.plot(seq(from=-20, to=30, by=2),seq(from=28, to=60, by=2),var.slice, 
           xlab = "Longitude", ylab = "Latitude", legend.shrink = 0.9, 
           legend.width = 1.2, legend.mar = NULL, main =nom_var,
           graphics.reset = FALSE, horizontal = FALSE, bigplot = NULL
           , smallplot = NULL, legend.only = FALSE, col = tim.colors(nlevel),
           lab.breaks=NULL, axis.args=NULL)


map(add=TRUE, fill=TRUE, col = alpha("grey", 0.5))
abline(v=seq(-21, 31, 2), lty=2)
abline(h=seq(27, 61, 2), lty=2)

Upvotes: 2

Views: 2446

Answers (1)

SESman
SESman

Reputation: 238

I've never used these functions but it seems rigorous to do so. Raster package offers some powerfull tools to reshape the data. When I worked (a few) with ncdf files I went directly to image() or image.plot{fields} to plot the map. With image.plot() you can give x and y values correponding to rows and columns to get good axis tick labels. I think it applies expand.grid itself. Once the ticks are good you can benefit from the tools of maps package easily.

Upvotes: 1

Related Questions