Here we attempt to provide useful code to generate figures from WRF outputs based on known galleries. For instance, NCL and WRF-Python provides extensive examples for plotting WRF outputs. Therefore, we aim to replicate some of these. Our approach to read wrfout files is based on eixport which relies r packages with GDAL bindings such as raster and stars. We do not try to provide a full gallery, instead, some basics and necessary plots to inspire other R used and receive more examples so share with the community.


  • eixport read and manipulate wrf files.
  • raster for gridded and raster data.
  • stars for gridded and raster data.
  • cptcity colour palettes.
  • sf for spatial vector data.
#> Loading required package: sp
#> Loading required package: abind
#> Loading required package: sf
#> Linking to GEOS 3.9.0, GDAL 3.2.2, PROJ 7.2.1

Based on NCL:

wrfo <- "/media/sergio/ext41/wrfo_sebrazil/wrfout_d02_2014-11-03_000000"

Reading HGT from wrfout

HGT <- wrf_get(wrfo,  "HGT", as_raster = T)
#> Warning in wrf_get(wrfo, "HGT", as_raster = T): The option as_raster will be
#> deprecated, see eixport::wrf_raster()
HGT <- HGT[[1]] # by default one variable with each time, so we select one

Adding coastlines and cropping for our study area

cl <- ne_countries(scale = "small", returnclass = "sf")
#> Warning in fun(libname, pkgname): rgeos: versions of GEOS runtime 3.9.0-CAPI-1.16.2
#> and GEOS at installation 3.8.0-CAPI-1.13.1differ
cl <- st_transform(cl, 31983)
HGT <- st_transform(st_as_stars(HGT), 31983)
cl <- st_cast(st_crop(cl, HGT), "LINESTRING")
#> Warning: attribute variables are assumed to be spatially constant throughout all
#> geometries
#> Warning in st_cast.sf(st_crop(cl, HGT), "LINESTRING"): repeating attributes for
#> all sub-geometries for which they may not be constant

Find colour palette for elevation

#> [1] "gmt_GMT_elevation" "grass_elevation"
#spplot(HGT, main = "HGT using spplot", scales=list(draw = TRUE),
#       col.regions = cpt("grass_elevation"),
#       sp.layout = list("sp.lines", as_Spatial(cl), col = "black"))
cols <- classInt::classIntervals(HGT$HGT_2014.11.03_00.00.00, n = 100, style = "pretty")
plot(HGT, axes = T, main = "HGT using stars", col = cpt("grass_elevation", n = length(cols$brks) - 1), breaks = cols$brks, reset = F)
plot(cl$geometry, add= T, col = "black")