SubsetNOMADS: Subset model data output

In rNOMADS: An R Interface to the NOAA Operational Model Archive and Distribution System

Subset model data output

Description

Subset data structure returned by ReadGrib or DODSGrab by variables, levels, etc.

Usage

SubsetNOMADS(model.data, levels = NULL, variables = NULL, lon = NULL, 
lat = NULL, ensembles = NULL, forecast.date = NULL, model.run.date = NULL)

Arguments

model.data	Data structure from ReadGrib or DODSGrab
levels	Vector of levels to keep
variables	Vector of variables to keep
lon	Vector of longitudes of model nodes to keep
lat	Vector of latitudes of model nodes to keep
ensembles	Vector of ensemble runs to keep
forecast.date	Vector of forecast dates to keep
model.run.date	Vector of model run dates to keep

Value

model.data.sub

A subset of model.data.

Note

Multiple elements in each argument vector are obviously OR (i.e. variables are “tmpprs” OR “hgtprs”) but multiple subset vectors are AND. Thus it is simple to construct a model.data.sub with variables: tmpprs and hgptprs only from ensemble runs 3 and 4, for example.

Author(s)

Daniel C. Bowman [e-mail]

See Also

ReadGrib, DODSGrab, ModelGrid, BuildProfile

Examples

## Not run: 
#Plot winds from 20 GENS model runs

#Get the latest ensemble model run
model.urls <- GetDODSDates("gens")
latest.model <- tail(model.urls$url, 1)
model.runs <- GetDODSModelRuns(latest.model)
model.run <- tail(model.runs$model.run[grepl("all", 
    model.runs$model.run)], 1)

#Define region of interest: Chapel Hill, NC
lon <- -79.052104
lat <- 35.907553

lons <- seq(0, 359, by = 1)
lats <- seq(-90, 90, by = 1)

lon.diff <- abs(lon + 360 - lons)
lat.diff <- abs(lat - lats)

model.lon.ind <- which(lon.diff == min(lon.diff)) - 1
model.lat.ind <- which(lat.diff == min(lat.diff)) - 1

#Set up call to NOMADS
time <- c(0, 0)      #Analysis(?) model only
node.lon  <- c(model.lon.ind - 2, model.lon.ind + 2)  #Longitude grid
node.lat  <- c(model.lat.ind - 2, model.lat.ind + 2)   #Latitude grid
variables <- c("ugrdprs", "vgrdprs", "hgtprs") #Wind speeds, and geopotential height
levels    <- c(0, 25) #All available levels
ensembles <- c(0, 20)  #All available ensembles

model.data <- DODSGrab(latest.model, model.run, variables, time, 
    node.lon, node.lat, levels = levels, ensembles = ensembles)

#Plot winds
zonal.wind <- NULL
merid.wind <- NULL
height     <- NULL

for(k in ((ensembles[1]:ensembles[2] + 1))) {
    model.data.sub <- SubsetNOMADS(model.data, ensembles = c(k), 
    variables = c("hgtprs", "ugrdprs", "vgrdprs"))
    profile <- BuildProfile(model.data.sub, lon + 360, lat)
    hgt     <- profile[[1]]$profile.data[, which(profile[[1]]$variables == "hgtprs"),]
    ugrd    <- profile[[1]]$profile.data[, which(profile[[1]]$variables == "ugrdprs"),]
    vgrd    <- profile[[1]]$profile.data[, which(profile[[1]]$variables == "vgrdprs"),]

   synth.hgt <- seq(min(hgt),
       max(hgt), length.out = 1000)
   ugrd.spline <- splinefun(hgt, ugrd, method = "natural")
   vgrd.spline <- splinefun(hgt, vgrd, method = "natural")
   zonal.wind[[k]] <- ugrd.spline(synth.hgt)
   merid.wind[[k]] <- vgrd.spline(synth.hgt)
   height[[k]] <- synth.hgt
}

PlotWindProfile(zonal.wind, merid.wind, height, lines = TRUE,
    points = FALSE, elev.circles = c(0, 15000, 30000), elev.labels = c(0, 15, 30),
    radial.lines = seq(45, 360, by = 45), colorbar = TRUE, invert = FALSE,
    point.cex = 2, pch = 19, lty = 1, lwd = 1,
    height.range = c(0, 30000), colorbar.label = "Wind Speed (m/s)")

## End(Not run)

rNOMADS documentation built on May 29, 2024, 6:44 a.m.