Note
Go to the end to download the full example code.
UserScript and StatAnalysis: Compute Polar Cap Temperature and Polar Vortex U and Create Plots
model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar.py
Scientific Objective
Many common modes of variability in the troposphere have stratospheric teleconnection pathways. Thus, the stratosphere can be a source of predictability for surface weather. This use case calls functions in METcalcpy to create polar cap temperature and polar vortex wind. It then runs Stat-Analysis on the output zonal means and creates a contour plot of bias in lead time and pressure level.
Version Added
METplus version 6.0
Datasets
Forecast: GFS Forecast U and T at multiple pressure levels
Observation: ERA Reanlaysis U and T at multiple pressure levels
Climatology: None
Location: Data for this use case is not contained in the sample data tar files due to its size. Rather, it is stored as additional data in a separate tar file, named additional_data_UserScript_fcstGFS_obsERA_StratosphereQBO.tar.gz and can be downloaded at https://dtcenter.ucar.edu/dfiles/code/METplus/METplus_Data/v6.0/.
METplus Components
This use case runs UserScript twice and Stat-Analysis once. The UserScripts compute polar cap temperature and polar vortex wind and create plots. METcalcpy, METplotpy, and METdataio are needed for this use case. The METcalcpy scripts accessed include the following:
metcalcpy/pre_processing/directional_means.py
metcalcpy/util/write_mpr.py
The METplotpy scripts accessed include the following:
metplotpy/contributed/stratosphere_diagnostics/stratosphere_plots.py
The METdataio scripts accessed include the following:
METreadnc/util/read_netcdf.py
METplus Workflow
Beginning time (VALID_BEG): 2018-02-01
End time (VALID_END): 2018-02-28
Increment between beginning and end times (VALID_INCREMENT): 30 days
Sequence of forecast leads to process (LEAD_SEQ): 0 - 384 hours at 24 hour intervals
This use case loops over lead times for the first UserScript and Stat-Analysis, and the second UserScript runs once over the entire time period. The first UserScript runs polar_t_u_driver.py which computes polar cap temperature and polar vortex wind and outputs that data to .stat files in MET’s matched pair format. Then, Stat-Analysis is run to compute bias and RMSE on the polar cap temperature and polar vortex wind. Finally, the second UserScript runs the bias_rmse_plot_driver.py which creates plots of bias and RMSE over lead time and pressure level.
METcalcpy 3.0.0 or higher, METplotpy 3.0.0 or higher, and METdataio 2.1 or higher are needed for this use case
METplus Configuration
METplus first loads all of the configuration files found in parm/metplus_config, then it loads any configuration files passed to METplus via the command line, i.e. parm/use_cases/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar.conf
[config]
# Documentation for this use case can be found at
# https://metplus.readthedocs.io/en/latest/generated/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar.html
# For additional information, please see the METplus Users Guide.
# https://metplus.readthedocs.io/en/latest/Users_Guide
###
# Processes to run
# https://metplus.readthedocs.io/en/latest/Users_Guide/systemconfiguration.html#process-list
###
PROCESS_LIST = UserScript(means), StatAnalysis(sanal_cnt), UserScript(plots_t), UserScript(plots_u)
#PROCESS_LIST = UserScript(means)
SCRUB_STAGING_DIR = False
###
# Time Info
# LOOP_BY options are INIT, VALID, RETRO, and REALTIME
# If set to INIT or RETRO:
# INIT_TIME_FMT, INIT_BEG, INIT_END, and INIT_INCREMENT must also be set
# If set to VALID or REALTIME:
# VALID_TIME_FMT, VALID_BEG, VALID_END, and VALID_INCREMENT must also be set
# LEAD_SEQ is the list of forecast leads to process
# https://metplus.readthedocs.io/en/latest/Users_Guide/systemconfiguration.html#timing-control
###
LOOP_BY = VALID
VALID_TIME_FMT = %Y%m%d
VALID_BEG = 20180201
VALID_END = 20180228
VALID_INCREMENT = 30d
#LEAD_SEQ = begin_end_incr(0,240,3),begin_end_incr(252,384,12)
LEAD_SEQ = begin_end_incr(0,384,24)
#LEAD_SEQ = 0
LOOP_ORDER = processes
###
# UserScript Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#userscript
###
[means]
USER_SCRIPT_RUNTIME_FREQ = RUN_ONCE_PER_LEAD
# Template of filenames to input to the user-script
USER_SCRIPT_INPUT_TEMPLATE = {INPUT_BASE}/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/ERA/ERA_{valid?fmt=%Y}_{valid?fmt=%m}.nc, {INPUT_BASE}/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/GFS/GFS_{valid?fmt=%Y}_{valid?fmt=%m}_{lead?fmt=%HHH}h.nc
# Name of the file containing the listing of input files
# The options are OBS_INPUT for observations or FCST_INPUT for forecast
# Or, set OBS_INPUT, FCST_INPUT if doing both and make sure the USER_SCRIPT_INPUT_TEMPLATE is ordered:
# observation_template, forecast_template
USER_SCRIPT_INPUT_TEMPLATE_LABELS = OBS_INPUT, FCST_INPUT
USER_SCRIPT_COMMAND = {PARM_BASE}/use_cases/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/polar_t_u_driver.py {lead?fmt=%HHH}
[user_env_vars]
MODEL_NAME = GFS
FCST_T_VAR = T
FCST_U_VAR = u
FCST_TIME_VAR = time
FCST_LAT_DIM = lat
FCST_LON_DIM = lon
FCST_PRES_DIM = pres
OBS_T_VAR = T
OBS_U_VAR = u
OBS_TIME_VAR = time
OBS_LAT_DIM = lat
OBS_LON_DIM = lon
OBS_PRES_DIM = pres
OUTPUT_DIR = {OUTPUT_BASE}/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar
PLOT_OUTPUT_DIR = {OUTPUT_DIR}/plots
PLOT_T_BIAS_LEVELS = -7,-6,-5,-4,-3,-2,-1,0,1,2,3,4,5,6,7,8,9
PLOT_T_BIAS_TITLE = GFS vs ERA Polar Cap Temperature Bias (ME) 02/2018
PLOT_T_RMSE_LEVELS = 0,1,2,3,4,5,6,7,8,9,10
PLOT_T_RMSE_TITLE = GFS vs ERA Polar Cap Temperature RMSE 02/2018
PLOT_T_BIAS_OUTPUT_FILE = ME_2018_02_polar_cap_T.png
PLOT_T_RMSE_OUTPUT_FILE = RMSE_2018_02_polar_cap_T.png
PLOT_U_BIAS_LEVELS = -4,-3,-2,-1,0,1,2,3,4,5,6,7,8,9,10,11,12
PLOT_U_BIAS_TITLE = GFS vs ERA Polar Vortex U Bias (ME) 02/2018
PLOT_U_RMSE_LEVELS = 0,2,4,6,8,10,12,14,16,18,20,22,24
PLOT_U_RMSE_TITLE = GFS vs ERA Polar Vortex U RMSE 02/2018
PLOT_U_BIAS_OUTPUT_FILE = ME_2018_02_polar_vortex_U.png
PLOT_U_RMSE_OUTPUT_FILE = RMSE_2018_02_polar_vortex_U.png
[sanal_cnt]
MODEL1 = GFS
MODEL1_OBTYPE = ADPUPA
STAT_ANALYSIS_CONFIG_FILE = {PARM_BASE}/met_config/STATAnalysisConfig_wrapped
STAT_ANALYSIS_JOB1 = -job aggregate_stat -line_type MPR -out_line_type CNT -fcst_var PolarCapT -by FCST_LEV,FCST_LEAD -out_stat [out_stat_file]_PolarCapT_CNT.stat
STAT_ANALYSIS_JOB2 = -job aggregate_stat -line_type MPR -out_line_type CNT -fcst_var PolarVortexU -by FCST_LEV,FCST_LEAD -out_stat [out_stat_file]_PolarVortexU_CNT.stat
MODEL_LIST = {MODEL1}
FCST_LEAD_LIST = 00, 03, 06, 09, 12, 15, 18, 21, 24, 27, 30, 33, 36, 39, 42, 45, 48, 51, 54, 57, 60, 63, 66, 69, 72, 75, 78, 81, 84, 87, 90, 93, 96, 99, 102, 105, 108, 111, 114, 117, 120, 123, 126, 129, 132, 135, 138, 141, 144, 147, 150, 153, 156, 159, 162, 165, 168, 171, 174, 177, 180, 183, 186, 189, 192, 195, 198, 201, 204, 207, 210, 213, 216, 219, 222, 225, 228, 231, 234, 237, 240, 252, 264, 276, 288, 300, 312, 324, 336, 348, 360, 372, 384
GROUP_LIST_ITEMS = MODEL_LIST
LOOP_LIST_ITEMS = FCST_LEAD_LIST
MODEL1_STAT_ANALYSIS_LOOKIN_DIR = {OUTPUT_BASE}/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/mpr
STAT_ANALYSIS_OUTPUT_DIR = {OUTPUT_BASE}/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/StatAnalysis
MODEL1_STAT_ANALYSIS_OUT_STAT_TEMPLATE = {model?fmt=%s}_ERA_{obs_valid_beg?fmt=%Y%m%d}_{obs_valid_end?fmt=%Y%m%d}_{lead?fmt=%H%M%S}L
[plots_t]
USER_SCRIPT_RUNTIME_FREQ = RUN_ONCE
USER_SCRIPT_INPUT_TEMPLATE = {OUTPUT_BASE}/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/StatAnalysis/GFS_ERA_{valid?fmt=%Y%m%d}_*_{lead?fmt=%H%M%S}L_PolarCapT_CNT.stat
# Name of the file containing the listing of input files
# The options are OBS_INPUT for observations or FCST_INPUT for forecast
# Or, set OBS_INPUT, FCST_INPUT if doing both and make sure the USER_SCRIPT_INPUT_TEMPLATE is ordered:
# observation_template, forecast_template
USER_SCRIPT_INPUT_TEMPLATE_LABELS = STAT_INPUT
USER_SCRIPT_COMMAND = {PARM_BASE}/use_cases/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/bias_rmse_plot_driver.py T
[plots_u]
USER_SCRIPT_RUNTIME_FREQ = RUN_ONCE
USER_SCRIPT_INPUT_TEMPLATE = {OUTPUT_BASE}/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/StatAnalysis/GFS_ERA_{valid?fmt=%Y%m%d}_*_{lead?fmt=%H%M%S}L_PolarVortexU_CNT.stat
# Name of the file containing the listing of input files
# The options are OBS_INPUT for observations or FCST_INPUT for forecast
# Or, set OBS_INPUT, FCST_INPUT if doing both and make sure the USER_SCRIPT_INPUT_TEMPLATE is ordered:
# observation_template, forecast_template
USER_SCRIPT_INPUT_TEMPLATE_LABELS = STAT_INPUT
USER_SCRIPT_COMMAND = {PARM_BASE}/use_cases/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/bias_rmse_plot_driver.py U
MET Configuration
METplus sets environment variables based on user settings in the METplus configuration file. See How METplus controls MET config file settings for more details.
YOU SHOULD NOT SET ANY OF THESE ENVIRONMENT VARIABLES YOURSELF! THEY WILL BE OVERWRITTEN BY METPLUS WHEN IT CALLS THE MET TOOLS!
If there is a setting in the MET configuration file that is currently not supported by METplus you’d like to control, please refer to: Overriding Unsupported MET config file settings
STATAnalysisConfig_wrapped
////////////////////////////////////////////////////////////////////////////////
//
// STAT-Analysis configuration file.
//
// For additional information, see the MET_BASE/config/README file.
//
////////////////////////////////////////////////////////////////////////////////
//
// Filtering input STAT lines by the contents of each column
//
//model = [
${METPLUS_MODEL}
//desc = [
${METPLUS_DESC}
//fcst_lead = [
${METPLUS_FCST_LEAD}
//obs_lead = [
${METPLUS_OBS_LEAD}
//fcst_valid_beg =
${METPLUS_FCST_VALID_BEG}
//fcst_valid_end =
${METPLUS_FCST_VALID_END}
fcst_valid_inc = [];
fcst_valid_exc = [];
//fcst_valid_hour = [
${METPLUS_FCST_VALID_HOUR}
//obs_valid_beg =
${METPLUS_OBS_VALID_BEG}
//obs_valid_end =
${METPLUS_OBS_VALID_END}
obs_valid_inc = [];
obs_valid_exc = [];
//obs_valid_hour = [
${METPLUS_OBS_VALID_HOUR}
//fcst_init_beg =
${METPLUS_FCST_INIT_BEG}
//fcst_init_end =
${METPLUS_FCST_INIT_END}
fcst_init_inc = [];
fcst_init_exc = [];
//fcst_init_hour = [
${METPLUS_FCST_INIT_HOUR}
//obs_init_beg =
${METPLUS_OBS_INIT_BEG}
//obs_init_end =
${METPLUS_OBS_INIT_END}
obs_init_inc = [];
obs_init_exc = [];
//obs_init_hour = [
${METPLUS_OBS_INIT_HOUR}
//fcst_var = [
${METPLUS_FCST_VAR}
//obs_var = [
${METPLUS_OBS_VAR}
//fcst_units = [
${METPLUS_FCST_UNITS}
//obs_units = [
${METPLUS_OBS_UNITS}
//fcst_lev = [
${METPLUS_FCST_LEVEL}
//obs_lev = [
${METPLUS_OBS_LEVEL}
//obtype = [
${METPLUS_OBTYPE}
//vx_mask = [
${METPLUS_VX_MASK}
//interp_mthd = [
${METPLUS_INTERP_MTHD}
//interp_pnts = [
${METPLUS_INTERP_PNTS}
//fcst_thresh = [
${METPLUS_FCST_THRESH}
//obs_thresh = [
${METPLUS_OBS_THRESH}
//cov_thresh = [
${METPLUS_COV_THRESH}
//alpha = [
${METPLUS_ALPHA}
//line_type = [
${METPLUS_LINE_TYPE}
column = [];
weight = [];
////////////////////////////////////////////////////////////////////////////////
//
// Array of STAT-Analysis jobs to be performed on the filtered data
//
//jobs = [
${METPLUS_JOBS}
////////////////////////////////////////////////////////////////////////////////
//
// Confidence interval settings
//
out_alpha = 0.05;
boot = {
interval = PCTILE;
rep_prop = 1.0;
n_rep = 0;
rng = "mt19937";
seed = "";
}
////////////////////////////////////////////////////////////////////////////////
//
// WMO mean computation logic
//
wmo_sqrt_stats = [ "CNT:FSTDEV", "CNT:OSTDEV", "CNT:ESTDEV",
"CNT:RMSE", "CNT:RMSFA", "CNT:RMSOA",
"VCNT:FS_RMS", "VCNT:OS_RMS", "VCNT:RMSVE",
"VCNT:FSTDEV", "VCNT:OSTDEV" ];
wmo_fisher_stats = [ "CNT:PR_CORR", "CNT:SP_CORR",
"CNT:KT_CORR", "CNT:ANOM_CORR" ];
////////////////////////////////////////////////////////////////////////////////
//hss_ec_value =
${METPLUS_HSS_EC_VALUE}
rank_corr_flag = FALSE;
vif_flag = FALSE;
tmp_dir = "${MET_TMP_DIR}";
//version = "V10.0";
${METPLUS_MET_CONFIG_OVERRIDES}
Python Embedding
This use case does not use python embedding.
User Scripting
There are two Python scripts in the is use case. The first, polar_t_u_driver.py reads in netCDF files for the forecast and observations and computes zonal means for temperature and wind. Then, polar cap temperature is computed by taking the meridional mean between 60 and 90 latitude, while polar vortex u is computed by taking the meridional mean between 50 and 80 latitude. This data is written to output files in MET’s matched pair format.
The second python script, bias_rmse_plot_driver.py reads in the output of Stat-Analysis and creates plots of bias and RMSE for polar cap temperature and polar vortex u. Variables input to this script are given in the [user_env_vars] section of the configuration file.
parm/use_cases/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/polar_t_u_driver.py
#!/usr/bin/env python3
"""
Create Polar Cap Temperatures and Polar Vortex Winds
"""
import os
import sys
import datetime
import numpy as np
import metcalcpy.pre_processing.directional_means as directional_means
import METreadnc.util.read_netcdf as read_netcdf
from metcalcpy.util.write_mpr import write_mpr_file
def main():
"""
Read arguments
"""
leadvar = sys.argv[1]
"""
Read METplus filename lists
"""
obs_filetxt = os.environ.get('METPLUS_FILELIST_OBS_INPUT','')
fcst_filetxt = os.environ.get('METPLUS_FILELIST_FCST_INPUT','')
with open(obs_filetxt) as ol:
obs_infiles = ol.read().splitlines()
# Remove the first line if it's there
if (obs_infiles[0] == 'file_list'):
obs_infiles = obs_infiles[1:]
with open(fcst_filetxt) as fl:
fcst_infiles = fl.read().splitlines()
# Remove the first line if it's there
if (fcst_infiles[0] == 'file_list'):
fcst_infiles = fcst_infiles[1:]
output_dir = os.environ['OUTPUT_DIR']
full_output_dir = os.path.join(output_dir,'mpr')
"""
Read variable/dimension names
"""
obs_tvar = os.environ.get('OBS_T_VAR','T')
obs_uvar = os.environ.get('OBS_U_VAR','u')
obs_latvar = os.environ.get('OBS_LAT_VAR','latitude')
obs_timevar = os.environ.get('OBS_TIME_VAR','time')
obs_latdim = os.environ.get('OBS_LAT_DIM','lat')
obs_londim = os.environ.get('OBS_LON_DIM','lon')
obs_presdim = os.environ.get('OBS_PRES_DIM','pres')
fcst_tvar = os.environ.get('FCST_T_VAR','T')
fcst_uvar = os.environ.get('FCST_U_VAR','u')
fcst_latvar = os.environ.get('FCST_LAT_VAR','latitude')
fcst_timevar = os.environ.get('FCST_TIME_VAR','time')
fcst_latdim = os.environ.get('FCST_LAT_DIM','lat')
fcst_londim = os.environ.get('FCST_LON_DIM','lon')
fcst_presdim = os.environ.get('FCST_PRES_DIM','pres')
"""
Read dataset
"""
file_reader = read_netcdf.ReadNetCDF()
dsO = file_reader.read_into_xarray(obs_infiles)[0]
dsO = dsO.rename({obs_latdim:'latitude',
obs_londim:'longitude',
obs_presdim:'pres'})
file_reader2 = read_netcdf.ReadNetCDF()
dsF = file_reader2.read_into_xarray(fcst_infiles)[0]
dsF = dsF.rename({fcst_latdim:'latitude',
fcst_londim:'longitude',
fcst_presdim:'pres'})
"""
Limit Dataset to 100 - 1 hPa
"""
dsO = dsO.sel(pres=slice(1,100))
dsF = dsF.sel(pres=slice(1,100))
"""
Assign Latitude Coordinate since it doesn't work
"""
dsO = dsO.assign_coords({obs_latvar:dsO[obs_latvar].values})
dsF = dsF.assign_coords({fcst_latvar:dsF[fcst_latvar].values})
"""
Create Polar Cap Temparatures for Forecast and Obs
"""
dsO = dsO.assign_coords({obs_latvar:dsO[obs_latvar].values})
TzmO = directional_means.zonal_mean(dsO[obs_tvar])
TO_6090 = directional_means.meridional_mean(TzmO, 60, 90)
TzmF = directional_means.zonal_mean(dsF[fcst_tvar])
TF_6090 = directional_means.meridional_mean(TzmF, 60, 90)
"""
Create Polar Vortex Winds
"""
UzmO = directional_means.zonal_mean(dsO[obs_uvar])
UO_6090 = directional_means.meridional_mean(UzmO, 50, 80)
UzmF = directional_means.zonal_mean(dsF[fcst_uvar])
UF_6090 = directional_means.meridional_mean(UzmF, 50, 80)
"""
Add P to the levels since they are pressure levels
"""
obs_lvls = ['P'+str(int(op)) for op in dsO['pres'].values]
obs_lvls2 = [str(int(op)) for op in dsO['pres'].values]
fcst_lvls = ['P'+str(int(fp)) for fp in dsF['pres'].values]
"""
Write output MPR files
"""
dlength1 = len(TO_6090[0,:])
dlength = dlength1*2
modname = os.environ.get('MODEL_NAME','GFS')
maskname = os.environ.get('MASK_NAME','FULL')
datetimeindex = dsF.indexes[fcst_timevar]
for i in range(len(datetimeindex)):
valid_str = datetimeindex[i].strftime('%Y%m%d_%H%M%S')
leadstr = str(int(leadvar)).zfill(2)+'0000'
outobs = np.concatenate((TO_6090[i,:].values,UO_6090[i,:].values))
outfcst = np.concatenate((TF_6090[i,:].values,UF_6090[i,:].values))
write_mpr_file(outfcst,outobs,[0.0]*dlength,[0.0]*dlength,[leadstr]*dlength,[valid_str]*dlength,
['000000']*dlength,[valid_str]*dlength,modname,['NA']*dlength,['PolarCapT']*dlength1+['PolarVortexU']*dlength1,
['K']*dlength1+['m/s']*dlength1,fcst_lvls*2,['PolarCapT']*dlength1+['PolarVortexU']*dlength1,
['K']*dlength1+['m/s']*dlength1,obs_lvls*2,maskname,obs_lvls2*2,full_output_dir,'polar_cap_T_stat_'+modname)
if __name__ == '__main__':
main()
parm/use_cases/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar/bias_rmse_plot_driver.py
#!/usr/bin/env python3
"""
Creates Polar Cap Bias and RMSE plots
"""
import os
import sys
import pandas as pd
import numpy as np
import METreadnc.util.read_netcdf as read_netcdf
from metplotpy.contributed.stratosphere_diagnostics.stratosphere_plots import plot_polar_bias,plot_polar_rmse
def main():
plvar = sys.argv[1]
"""
Read METplus environment variables
"""
plot_bias_levels_str = os.environ['PLOT_'+plvar+'_BIAS_LEVELS'].split(',')
plot_bias_levels = [float(pp) for pp in plot_bias_levels_str]
plot_bias_title = os.environ['PLOT_'+plvar+'_BIAS_TITLE']
plot_rmse_levels_str = os.environ['PLOT_'+plvar+'_RMSE_LEVELS'].split(',')
plot_rmse_levels = [float(pp) for pp in plot_rmse_levels_str]
plot_rmse_title = os.environ['PLOT_'+plvar+'_RMSE_TITLE']
output_dir = os.environ['PLOT_OUTPUT_DIR']
bias_output_file = os.environ.get('PLOT_'+plvar+'_BIAS_OUTPUT_FILE','bias_plot.png')
rmse_output_file = os.environ.get('PLOT_'+plvar+'_RMSE_OUTPUT_FILE','rmse_plot.png')
plot_bias_output_file = os.path.join(output_dir,bias_output_file)
plot_rmse_output_file = os.path.join(output_dir,rmse_output_file)
"""
Make plot output directory if it doesn't exist
"""
if not os.path.exists(output_dir):
os.makedirs(output_dir)
"""
Read the list of files
"""
stat_filetxt = os.environ.get('METPLUS_FILELIST_STAT_INPUT','')
with open(stat_filetxt) as sf:
stat_infiles = sf.read().splitlines()
# Remove the first line if it's there
if (stat_infiles[0] == 'file_list'):
stat_infiles = stat_infiles[1:]
"""
Read the first file and set up arrays
"""
fleads = len(stat_infiles)
dfin = pd.DataFrame(pd.read_csv(stat_infiles[0],delim_whitespace=True,header=0))
dfin['FCST_LEV'] = dfin['FCST_LEV'].str.replace('P', '')
dfin['FCST_LEV'] = dfin['FCST_LEV'].astype('float64')
dfin = dfin.sort_values('FCST_LEV')
flvls = len(dfin)
plevels = np.empty([fleads,flvls],dtype=float)
pleads = np.empty([fleads,flvls],dtype=float)
prmse = np.empty([fleads,flvls],dtype=float)
pme = np.empty([fleads,flvls],dtype=float)
plevels[0,:] = dfin['FCST_LEV'].to_numpy()
pleads[0,:] = dfin['FCST_LEAD'].to_numpy()/10000
prmse[0,:] = dfin['RMSE'].astype('float64')
pme[0,:] = dfin['ME'].astype('float64')
"""
Read in the rest of the data
"""
for i in range(1,len(stat_infiles)):
df = pd.DataFrame(pd.read_csv(stat_infiles[i],delim_whitespace=True,header=0))
df['FCST_LEV'] = df['FCST_LEV'].str.replace('P', '')
df['FCST_LEV'] = df['FCST_LEV'].astype('float64')
dfnew = df.sort_values('FCST_LEV')
plevels[i,:] = dfnew['FCST_LEV'].to_numpy()
#pleads[i] = dfnew['FCST_LEAD'].to_numpy()[0]/10000
pleads[i,:] = dfnew['FCST_LEAD'].to_numpy()/10000
prmse[i,:] = dfnew['RMSE'].astype('float64')
pme[i,:] = dfnew['ME'].astype('float64')
"""
Create plots
"""
plot_polar_bias(pleads,plevels,pme,plot_bias_output_file,plot_bias_title,plot_bias_levels)
plot_polar_rmse(pleads,plevels,prmse,plot_rmse_output_file,plot_rmse_title,plot_rmse_levels)
if __name__ == '__main__':
main()
Running METplus
Pass the use case configuration file to the run_metplus.py script along with any user-specific system configuration files if desired:
run_metplus.py /path/to/METplus/parm/use_cases/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar.conf /path/to/user_system.conf
See Running METplus for more information.
Expected Output
A successful run will output the following both to the screen and to the logfile:
INFO: METplus has successfully finished running.
Refer to the value set for OUTPUT_BASE to find where the output data was generated. Output for this use case will be found in {OUTPUT_BASE}/model_applications/s2s_stratosphere/UserScript_fcstGFS_obsERA_StratospherePolar The output includes matched pair files, one for each lead time and valid time in mpr directory. The output files have the following format:
* polar_cap_T_stat_GFS_HH0000L_YYYYMMDD_000000V.stat
* polar_cap_T_stat_GFS_HHH0000L_YYYYMMDD_000000V.stat
Additionally, CNT statistics will be output to the StatAnalysis directory with one file for each lead time and variable. The files have the following format:
* GFS_ERA_20180201_20180228_HH0000L_PolarCapT_CNT.stat
* GFS_ERA_20180201_20180228_HH0000L_PolarVortexU_CNT.stat
Four plot are also output to the plots directory:
* ME_2018_02_polar_cap_T.png
* ME_2018_02_polar_vortex_U.png
* RMSE_2018_02_polar_cap_T.png
* RMSE_2018_02_polar_vortex_U.png
Keywords
Note
S2SAppUseCase
S2SStratosphereAppUseCase
UserScriptUseCase
StatAnalysisUseCase
METdataioUseCase
METcalcpyUseCase
METplotpyUseCase
Navigate to the METplus Quick Search for Use Cases page to discover other similar use cases.
sphinx_gallery_thumbnail_path = ‘_static/s2s_stratosphere-UserScript_fcstGFS_obsERA_StratospherePolar.png’