Weather regimes can be calculated using the k-means clustering algorithm, typically applied to the PCs of anomalies from selected variables, after smoothing and detrending.
The work flow will be something like:
Load()
.Ano()
.Smooth()
.Trend()
.WeatherRegime()
.PlotEquiMap()
.The loading, detrending and calculation of anomalies needs to be done before applying the WeatherRegime() function. The calculation of the EOFs/PC (including the latitudinal area weighting) and the reconstruction of the field (Cluster centers * EOFs) is done within the WeatherRegime() function.
Below is a simple example of calculating the weather regimes for the North Atlantic from JRA55 SLP for January (without smoothing/detrending).
```r
library(s2dverification)
library(reshape2)
JRA55 <- '/esnas/recon/jma/jra55/$STORE_FREQ$_mean/$VAR_NAME$_f6h/$VAR_NAME$_$YEAR$$MONTH$.nc'
psl <- "psl"
year.start <- 1981
year.end <- 2016
#Euro-Atlantic region
lat.max <- 81
lat.min <- 27
lon.max <- 45
lon.min <- 274.5
firstmonth <- 1
leadmonth <- 0
psl <- Load(var = psl, exp = NULL, obs = list(list(path=JRA55)), paste0(year.start:year.end,'0101'), storefreq = 'daily', leadtimemax = 366, output = 'lonlat', latmin = lat.min, latmax = lat.max, lonmin = lon.min, lonmax = lon.max, nprocs=8)
lon <- psl$lon
lat <- psl$lat
psl <- psl$obs
#psl <- Trend(psl, posTR = 3)$detrended
psl <- psl[1, 1, , , , ] # Drop the dataset and members dimensions
psl <- Subset(psl, 2, 1 : 31) # Select the month of January
psl_ano <- psl - InsertDim(Mean1Dim(psl, 1), 1, 36)
names(dim(psl_ano)) <- c("ftime", "sdate", "lat", "lon")
WR <- WeatherRegime(psl_ano, lon = lon, lat = lat, ncenters = 4)
PlotEquiMap(WR$field[1,,], lon = lon, lat = lat) #Plot the 1st regime
```