getClassification: Getting the cost733 WLK Classification For a Given Date
In retostauffer/Rcost: Cost 733 WLK Classification in R

Description Usage Arguments Details Value Author(s) References Examples

This function was the reason to set up this small R package. This package performs the cost733 WLK classification on gridde data sets such as the ECMWF INTERIM or ECMWF HIRES. The classification is performed based on the fortran routine by Andreas Philipp, UNI Augsburg.

http://cost733class.geo.uni-augsburg.de/moin/cost733wiki/data/pages_bak/WLKC733/attachments/cost733_wlk_1.07.f90

The classification requires the variables “tcw” (total column water), “u700” and “v700” (zonal and meridional wind component in 700 hectopascal), as far as “z500” and “z925”. These variables have to be provided by the nc files given as input argument. See details on how the netcdf files should look like.

The defaults of this method use the specificationas as proposed by Thomas Krennert (ZAMG Vienna). The method here is not 100 percent generic as it e.g, expects a specific file format (and file names) was I have only planned to use the script on ECMWF 00UTC deterministic forecasts at the moment.

## Basic usage
getClassification(nc, date, steps, subset, pwclim, weights,
                  nsector=8, maindirthreshold=0.4,
                  zamg=FALSE, verbose=0)

`nc`	`list` containing one or more open netcdf connection opened by `ncdf4::nc_open`.
`date`	`Date`, `POSIXt`, or `character` object which can be converted into `POSIXt`. Date/time when the forecast was initialized. Used to find the NetCDF files with a specific (fixed) name at the moment.
`steps`	Optional, if not set the classification will be performed on all steps provided by the NetCDF files.
`subset`	Optional, default `NULL`. See `getdata` for more details.
`pwclim`	Optional. If not set the dry/wet classification cannot be performed. If set it has to be a `numeric` vector of length `366` containing the climatological “precip. water content” for each day of the year.
`weights`	Optional, if not set the ZAMG config with `c(15,2,1)` will be used for the weighting. If set it has to be a `numeric` vector of length 3 (for core/mid/margin weight).
`nsector`	Integer, default `8`. Number of wind sectors for the classification.
`maindirthreshold`	`numeric`, default is `0.4`. Used for wind sector classification. If less than `maindirthreshold` can be classified for one main wind direction the wind sector classification will be “variable” (`sector=0`).
`zamg`	Logical flag, default `FALSE`. If set to `TRUE` the main wind direction will be computed based on the classification by “Trendanalyse Klien Endbericht” (as suggested by Thomas Kennert, ZAMG).
`verbose`	Integer. Logical `TRUE`/`FALSE` will be converted to `0`/`1`. However, verbose levels up to `4` (highly verbose) can be set.

For subset please see details in getdata.

The netcdf files have to provide the variables u and v (wind components) on 700 hectopascal, z (geopotential height) on 500 and 925 hectopascal, and tcw (total column water). To be able to extract the necessary information the netcdf file(s) have to provide the following information: - dimension longitude (named longitude) - dimension latitude (named latitude) - dimension time (named time; units something like “hours since 1900-01-01 00:00” or similar) - dimension level (named level) which is used when searching for level variables I typically use two different files, a surface level netcdf (dimensions: longitude, latitude, and time; containing at least tcw), and a pressure level file (dimensions: longitude, latitude, time, and level; containing at least u700, v700, z500, z925) and use both as input arguments. E.g., - getClassification( list(nc_open('<sfcfile>'),nc_open('<plfile>')), ... )

Returns a zoo object with the following items: - date: date/time when the classification is valid - yday: day of the year - sector: classifier: wind sector - cyclonic925/cyclonic500: classifier: wheter it is 0=anticyclonic or 1=cyclonic in 925/500 hPa - wet: classifier: -99=missing (due to lack of pwclim), 0=dry, and 1=wet. - c925/c500: curvature values used to create the classification. - pwi: computed pwi values which will be used for classification

Reto Stauffer

Philipp A., C. Beck, R. Huth and J. Jacobeit (2014): Development and comparison of circulation type classifications using the COST 733 dataset and software. International Journal of Climatology. DOI: 10.1002/joc.3920

# Open netcdf connections
require("ncdf4")
files <- demofiles()
nc    <- list(nc_open(files$sfc), nc_open(files$pl)) # Open NetCDF Files

# Extract initial date (hours since 1900-01-01 00:00:0.0)
# Take first date only.
init  <- as.POSIXct( ncvar_get(nc[[1]],"time") * 3600, origin = "1900-01-01" )[1]

# Perform classification
x1 <- getClassification( nc, init)
plot(x1, main = "COST733 Classification Output (x1)")

# Use different areal weights and only 4 wind sectors
x2 <- getClassification( nc, init, weights = c(10,5,2), nsector = 4 )
plot(x2, main = "COST733 Classification Output (x2)")

# Provide pw climatology (in this case simply a constant of
# 24.844 for all 1:365 days of a year): returns additional
# binary classification "wet" (simply pwi > pwiclim).
x3 <- getClassification( nc, init, weights = c(10,5,2), nsector = 8,
                         pwclim = rep(24.844,365) )
plot(x3, main = "COST733 Classification Output (x3)")

# Close NetCDF Files
for ( n in nc ) nc_close(n)