Note
Go to the end to download the full example code
Blocking Calculation: GFS and ERA RegridDataPlane, PcpCombine, and Blocking python code
model_applications/ s2s_mid_lat/ UserScript_fcstGFS_obsERA_Blocking.py
Scientific Objective
To compute the Central Blocking Latitude, Instantaneousy blocked latitudes, Group Instantaneousy blocked latitudes, and the frequency of atmospheric blocking using the Pelly-Hoskins Method. After these are computed, contingency table statistics are computed on the Instantaneous blocked latitudes and blocks using stat_analysis.
Datasets
Forecast dataset: GFS Forecast 500 mb height.
Observation dataset: ERA Reanlaysis 500 mb height.
External Dependencies
You will need to use a version of Python 3.6+ that has the following packages installed:
* numpy
* netCDF4
* datetime
* bisect
* scipy
If the version of Python used to compile MET did not have these libraries at the time of compilation, you will need to add these packages or create a new Python environment with these packages.
If this is the case, you will need to set the MET_PYTHON_EXE environment variable to the path of the version of Python you want to use. If you want this version of Python to only apply to this use case, set it in the [user_env_vars] section of a METplus configuration file.:
[user_env_vars] MET_PYTHON_EXE = /path/to/python/with/required/packages/bin/python
METplus Components
This use case runs the blocking driver script which runs the steps the user lists in STEPS_OBS. The possible steps are regridding, time averaging, computing a running mean, computing anomalies, computing CBLs (CBL), plotting CBLs (PLOTCBL), computing IBLs (IBL), plotting IBL frequency (PLOTIBL), computing GIBLs (GIBL), computing blocks (CALCBLOCKS), plotting the blocking frequency (PLOTBLOCKS) and using stat_analysis to compute statistics on the IBL or blocking results. Regridding, time averaging, running means, anomaloies, and stat_analysis are set up in the UserScript .conf file and are formatted as follows: PROCESS_LIST = RegridDataPlane(regrid_fcst), RegridDataPlane(regrid_obs), PcpCombine(daily_mean_fcst), PcpCombine(daily_mean_obs), PcpCombine(running_mean_obs), PcpCombine(anomaly_obs), UserScript(create_cbl_filelist), UserScript(script_blocking), StatAnalysis(sanal_ibls), StatAnalysis(sanal_blocks)
The other steps are listed in the Blocking .conf file and are formatted as follows: FCST_STEPS = CBL+IBL+PLOTIBL+GILB+CALCBLOCKS+PLOTBLOCKS OBS_STEPS = CBL+PLOTCBL+IBL+PLOTIBL+GILB+CALCBLOCKS+PLOTBLOCKS
METplus Workflow
The blocking python code is run for each time for the forecast and observations data. This example loops by init time for the model pre-processing, and valid time for the other steps. This version is set to only process the blocking steps (CBL, PLOTCBL, IBL, PLOTIBL) and stat_analysis, omitting the regridding, time averaging, running mean, and anomaly pre-processing steps. However, the configurations for pre-processing are available for user reference.
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_mid_lat/UserScript_fcstGFS_obsERA_Blocking.py. The file UserScript_fcstGFS_obsERA_Blocking.conf runs the python program, and the variables for all steps of the Blocking calculation are given in the [user_env_vars] section of the .conf file.
[config]
# Documentation for this use case can be found at
# https://metplus.readthedocs.io/en/latest/generated/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking.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(create_cbl_filelist), UserScript(script_blocking), StatAnalysis(sanal_ibls), StatAnalysis(sanal_blocks)
# use this process list if pre-processing steps are needed
# PROCESS_LIST = RegridDataPlane(regrid_fcst), RegridDataPlane(regrid_obs), PcpCombine(daily_mean_fcst), PcpCombine(daily_mean_obs), PcpCombine(running_mean_obs), PcpCombine(anomaly_obs), UserScript(create_cbl_filelist), UserScript(script_blocking)
###
# 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%H
VALID_BEG = 2000120100
VALID_END = 2017022800
VALID_INCREMENT = 86400
LEAD_SEQ = 0
# Only Process DJF
SKIP_TIMES = "%m:begin_end_incr(3,11,1)", "%m%d:0229"
###
# File I/O
# https://metplus.readthedocs.io/en/latest/Users_Guide/systemconfiguration.html#directory-and-filename-template-info
###
OBS_ANOM_INPUT_DIR = {INPUT_BASE}/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/ERA/Anomaly
OBS_ANOM_INPUT_TEMPLATE = Z500_anomaly_{valid?fmt=%Y%m%d}_NH.nc
OBS_ANOM_OUTPUT_DIR = {OBS_ANOM_INPUT_DIR}
OBS_ANOM_OUTPUT_TEMPLATE = ERA_anom_files_lead{lead?fmt=%HHH}.txt
OBS_AVE_INPUT_DIR = {INPUT_BASE}/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/ERA/Daily
OBS_AVE_INPUT_TEMPLATE = Z500_daily_{valid?fmt=%Y%m%d}_NH.nc
OBS_AVE_OUTPUT_DIR = {OBS_AVE_INPUT_DIR}
OBS_AVE_OUTPUT_TEMPLATE = ERA_daily_files_lead{lead?fmt=%HHH}.txt
FCST_AVE_INPUT_DIR = {INPUT_BASE}/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/GFS/Daily
FCST_AVE_INPUT_TEMPLATE = Z500_daily_{init?fmt=%Y%m%d}_{lead?fmt=%HHH}_NH.nc
FCST_AVE_OUTPUT_DIR = {FCST_AVE_INPUT_DIR}
FCST_AVE_OUTPUT_TEMPLATE = GFS_daily_files_lead{lead?fmt=%HHH}.txt
###
# RegridDataPlane(regrid_fcst) Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#regriddataplane
###
# Forecast Regridding to 1 degree using regrid_data_plane
[regrid_fcst]
LOOP_BY = INIT
INIT_TIME_FMT = %Y%m%d%H
INIT_BEG = 2000120100
INIT_END = 2017022800
INIT_INCREMENT = 86400
LEAD_SEQ = 24
# REGRID_DATA_PLANE (Step 1)
FCST_REGRID_DATA_PLANE_RUN = True
FCST_DATA_PLANE_ONCE_PER_FIELD = False
FCST_REGRID_DATA_PLANE_VAR1_INPUT_FIELD_NAME = Z500
FCST_REGRID_DATA_PLANE_VAR1_INPUT_LEVEL = P500
FCST_REGRID_DATA_PLANE_VAR1_OUTPUT_FIELD_NAME = Z500
REGRID_DATA_PLANE_VERIF_GRID = latlon 360 90 89 0 -1.0 1.0
REGRID_DATA_PLANE_METHOD = BILIN
REGRID_DATA_PLANE_WIDTH = 2
FCST_REGRID_DATA_PLANE_INPUT_DIR = /gpfs/fs1/p/ral/jntp/GMTB/Phys_Test_FV3GFSv2/POST/suite1/
FCST_REGRID_DATA_PLANE_OUTPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/FV3GFS/Regrid
FCST_REGRID_DATA_PLANE_INPUT_TEMPLATE = {init?fmt=%Y%m%d%H}/gfs.t00z.pgrb2.0p25.f{lead?fmt=%HHH}
FCST_REGRID_DATA_PLANE_OUTPUT_TEMPLATE = {init?fmt=%Y%m%d%H}/Z500_3hourly_{init?fmt=%Y%m%d%H}_{lead?fmt=%HHH}_NH.nc
###
# RegridDataPlane(regrid_obs) Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#regriddataplane
###
# Observation Regridding to 1 degree using regrid_data_plane
[regrid_obs]
LOOP_BY = VALID
VALID_TIME_FMT = %Y%m%d%H
VALID_BEG = 1979120100
VALID_END = 2017022818
VALID_INCREMENT = 21600
LEAD_SEQ = 0
OBS_REGRID_DATA_PLANE_RUN = True
OBS_DATA_PLANE_ONCE_PER_FIELD = False
OBS_REGRID_DATA_PLANE_VAR1_INPUT_FIELD_NAME = Z
OBS_REGRID_DATA_PLANE_VAR1_INPUT_LEVEL = P500
OBS_REGRID_DATA_PLANE_VAR1_OUTPUT_FIELD_NAME = Z500
REGRID_DATA_PLANE_VERIF_GRID = latlon 360 90 89 0 -1.0 1.0
REGRID_DATA_PLANE_METHOD = BILIN
REGRID_DATA_PLANE_WIDTH = 2
OBS_REGRID_DATA_PLANE_INPUT_DIR = /gpfs/fs1/collections/rda/data/ds627.0/ei.oper.an.pl
OBS_REGRID_DATA_PLANE_OUTPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/ERA/Regrid
OBS_REGRID_DATA_PLANE_INPUT_TEMPLATE = {valid?fmt=%Y%m}/ei.oper.an.pl.regn128sc.{valid?fmt=%Y%m%d%H}
OBS_REGRID_DATA_PLANE_OUTPUT_TEMPLATE = {valid?fmt=%Y%m%d}/Z500_6hourly_{init?fmt=%Y%m%d%H}_NH.nc
###
# PCPCombine(daily_mean_obs) Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#pcpcombine
###
# Perform a sum over the 4 daily times that have been regridded using pcp_combine
# 00, 06, 12, 18 UTC
[daily_mean_obs]
LOOP_BY = VALID
VALID_TIME_FMT = %Y%m%d%H
VALID_BEG = 1979120118
VALID_END = 2017022818
VALID_INCREMENT = 86400
OBS_PCP_COMBINE_RUN = True
OBS_PCP_COMBINE_METHOD = DERIVE
OBS_PCP_COMBINE_STAT_LIST = MEAN
OBS_PCP_COMBINE_INPUT_ACCUMS = 6
OBS_PCP_COMBINE_INPUT_NAMES = Z500
OBS_PCP_COMBINE_INPUT_LEVELS = "(*,*)"
OBS_PCP_COMBINE_INPUT_OPTIONS = convert(x) = x / 9.81; set_attr_valid = "{valid?fmt=%Y%m%d_%H%M%S?shift=-64800}";
OBS_PCP_COMBINE_OUTPUT_NAME = Z500
OBS_PCP_COMBINE_OUTPUT_ACCUM = 24
OBS_PCP_COMBINE_DERIVE_LOOKBACK = 24
OBS_PCP_COMBINE_INPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/ERA/Regrid
OBS_PCP_COMBINE_OUTPUT_DIR = {OBS_AVE_INPUT_DIR}
# Input ERA Interim
OBS_PCP_COMBINE_INPUT_TEMPLATE = {valid?fmt=%Y%m%d}/Z500_6hourly_{valid?fmt=%Y%m%d%H}_NH.nc
OBS_PCP_COMBINE_OUTPUT_TEMPLATE = {OBS_AVE_INPUT_TEMPLATE}
###
# PCPCombine(running_mean_obs) Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#pcpcombine
###
# Perform a 5 day running mean on the data using pcp_combine
[running_mean_obs]
LOOP_BY = VALID
VALID_TIME_FMT = %Y%m%d%H
VALID_BEG = 1979120100
VALID_END = 2017022800
VALID_INCREMENT = 86400
# Add the first/last 2 days to the skip times to compute the running mean
SKIP_TIMES = "%m:begin_end_incr(3,11,1)", "%m%d:1201,1202,0227,0228,0229"
OBS_PCP_COMBINE_RUN = TRUE
OBS_PCP_COMBINE_METHOD = DERIVE
OBS_PCP_COMBINE_STAT_LIST = MEAN
OBS_PCP_COMBINE_INPUT_ACCUMS = 24
OBS_PCP_COMBINE_INPUT_NAMES = Z500
OBS_PCP_COMBINE_INPUT_LEVELS = "(*,*)"
OBS_PCP_COMBINE_INPUT_OPTIONS = set_attr_valid = "{valid?fmt=%Y%m%d_%H%M%S?shift=-172800}";
OBS_PCP_COMBINE_OUTPUT_NAME = Z500
# Running mean is 5 days
OBS_PCP_COMBINE_OUTPUT_ACCUM = 120
OBS_PCP_COMBINE_DERIVE_LOOKBACK = 120
OBS_PCP_COMBINE_INPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/ERA/Daily
OBS_PCP_COMBINE_OUTPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/ERA/Rmean5d
OBS_PCP_COMBINE_INPUT_TEMPLATE = Z500_daily_{valid?fmt=%Y%m%d}_NH.nc
OBS_PCP_COMBINE_OUTPUT_TEMPLATE = Z500_5daymean_{valid?fmt=%Y%m%d?shift=-172800}_NH.nc
###
# PCPCombine(anomaly_obs) Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#pcpcombine
###
# Compute anomalies using the daily means and 5 day running mean using pcp_combine
[anomaly_obs]
LOOP_BY = VALID
VALID_TIME_FMT = %Y%m%d%H
VALID_BEG = 1979120100
VALID_END = 2017022800
VALID_INCREMENT = 86400
LEAD_SEQ = 0
# Add the first/last 2 days to the skip times to compute the running mean
SKIP_TIMES = "%m:begin_end_incr(3,11,1)", "%m%d:1201,1202,0227,0228,0229"
OBS_PCP_COMBINE_RUN = True
OBS_PCP_COMBINE_METHOD = USER_DEFINED
OBS_PCP_COMBINE_COMMAND = -subtract {OBS_PCP_COMBINE_INPUT_DIR}/Daily/Z500_daily_{valid?fmt=%Y%m%d}_NH.nc {OBS_PCP_COMBINE_INPUT_DIR}/Rmean5d/Z500_5daymean_{valid?fmt=%Y%m%d}_NH.nc -field 'name="Z500"; level="(*,*)";'
OBS_PCP_COMBINE_INPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/ERA
OBS_PCP_COMBINE_OUTPUT_DIR = {OBS_ANOM_INPUT_DIR}
OBS_PCP_COMBINE_INPUT_TEMPLATE = Z500_daily_{valid?fmt=%Y%m%d}_NH.nc
OBS_PCP_COMBINE_OUTPUT_TEMPLATE = {OBS_ANOM_INPUT_TEMPLATE}
###
# UserScript(create_cbl_filelist) Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#userscript
###
# This is run separately since it has different start/end times
[create_cbl_filelist]
# Skip the days on the edges that are not available due to the running mean
SKIP_TIMES = "%m:begin_end_incr(3,11,1)", "%m%d:1201,0229"
# Find the files for each lead time
USER_SCRIPT_RUNTIME_FREQ = RUN_ONCE_PER_LEAD
# Valid Begin and End Times for the CBL File Climatology
VALID_BEG = 1979120100
VALID_END = 2017022800
VALID_INCREMENT = 86400
LEAD_SEQ = 0
USER_SCRIPT_INPUT_TEMPLATE = {INPUT_BASE}/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/ERA/Anomaly/Z500_anomaly_{valid?fmt=%Y%m%d}_NH.nc
# Name of the file containing the listing of input files
USER_SCRIPT_INPUT_TEMPLATE_LABELS = OBS_CBL_INPUT
# Placeholder command just to build the file list
# This just states that it's building the file list
USER_SCRIPT_COMMAND = echo Populated file list for CBL Input
[user_env_vars]
# Obs and/or Forecast
FCST_STEPS = CBL+IBL+PLOTIBL+GIBL+CALCBLOCKS+PLOTBLOCKS
OBS_STEPS = CBL+PLOTCBL+IBL+PLOTIBL+GIBL+CALCBLOCKS+PLOTBLOCKS
# Number of Seasons and Days per season that should be available
# The code will fill missing data, but requires the same number of days per
# season for each year. You may need to omit leap days if February is part of
# the processing
CBL_NUM_SEASONS = 38
IBL_NUM_SEASONS = 17
DAYS_PER_SEASON = 89
# Use the obs climatology for the calculation of CBL data because the forecast
# does not have a long enough data history. Set to False if not wanting to
# use the obs
USE_CBL_OBS = True
# Variable Name for the Z500 anomaly data to read in to the blocking python code
OBS_BLOCKING_ANOMALY_VAR = Z500_ANA
# Variable for the Z500 data
FCST_BLOCKING_VAR = Z500_P500
OBS_BLOCKING_VAR = Z500
# Number of model grid points used for a moving average
# Must be odd
FCST_SMOOTHING_PTS = 9
OBS_SMOOTHING_PTS = {FCST_SMOOTHING_PTS}
# Lat Delta, to allow for offset from the Central Blocking Latitude
FCST_LAT_DELTA = -5,0,5
OBS_LAT_DELTA = {FCST_LAT_DELTA}
# Meridional Extent of blocks (NORTH_SOUTH_LIMITS/2)
FCST_NORTH_SOUTH_LIMITS = 30
OBS_NORTH_SOUTH_LIMITS = {FCST_NORTH_SOUTH_LIMITS}
# Maximum number of grid points between IBLs for everything in between to be included as an IBL
FCST_IBL_DIST = 7
OBS_IBL_DIST = {FCST_IBL_DIST}
# Number of grid points in and IBL to make a GIBL
FCST_IBL_IN_GIBL = 15
OBS_IBL_IN_GIBL = {FCST_IBL_IN_GIBL}
# Number of grid points that must overlap across days for a GIBL
FCST_GIBL_OVERLAP = 10
OBS_GIBL_OVERLAP = {FCST_GIBL_OVERLAP}
# Time duration in days needed for a block
FCST_BLOCK_TIME = 5
OBS_BLOCK_TIME = {FCST_BLOCK_TIME}
# Number of grid points a block must travel to terminate
FCST_BLOCK_TRAVEL = 45
OBS_BLOCK_TRAVEL = {FCST_BLOCK_TRAVEL}
# Method to compute blocking. Currently, the only option is 'PH' for the
# Pelly-Hoskins Method
FCST_BLOCK_METHOD = PH
OBS_BLOCK_METHOD = {FCST_BLOCK_METHOD}
# Location of output MPR files
BLOCKING_MPR_OUTPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/mpr
# Plots Output Dir
BLOCKING_PLOT_OUTPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/plots
#CBL plot title and output namename
OBS_CBL_PLOT_MTHSTR = DJF
OBS_CBL_PLOT_OUTPUT_NAME = ERA_CBL_avg
# IBL plot title and output name
IBL_PLOT_TITLE = DJF Instantaneous Blocked Longitude
IBL_PLOT_OUTPUT_NAME = FV3_ERA_IBL_Freq_DJF
# IBL plot legend for forecast and obs
IBL_PLOT_OBS_LABEL = ERA Reanalysis
IBL_PLOT_FCST_LABEL = GEFS
###
# UserScript(script_blocking) Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#userscript
###
# Run the Blocking Analysis Script
[script_blocking]
# Timing Information
LEAD_SEQ = 24
# Skip the days on the edges that are not available due to the running mean
SKIP_TIMES = "%m:begin_end_incr(3,11,1)", "%m%d:1201,0229"
# Run the user script once for each lead
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_mid_lat/UserScript_fcstGFS_obsERA_Blocking/ERA/Daily/Z500_daily_{valid?fmt=%Y%m%d}_NH.nc,{INPUT_BASE}/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/GFS/Daily/Z500_{init?fmt=%Y%m%d}_{lead?fmt=%HHH}_NH.nc
# Name of the file containing the listing of input files
# The options are OBS_CBL_INPUT, FCST_CBL_INPUT, OBS_IBL_INPUT, and FCST_IBL_INPUT
# *** Make sure the order is the same as the order of templates listed in USER_SCRIPT_INPUT_TEMPLATE
USER_SCRIPT_INPUT_TEMPLATE_LABELS = OBS_IBL_INPUT, FCST_IBL_INPUT
# Command to run the user script with input configuration file
USER_SCRIPT_COMMAND = {METPLUS_BASE}/parm/use_cases/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/Blocking_driver.py
###
# StatAnalysis(sanal_ibls) Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#statanalysis
###
# Stat Analysis for the IBLs
[sanal_ibls]
VALID_TIME_FMT = %Y%m%d
VALID_BEG = 20001201
VALID_END = 20170228
STAT_ANALYSIS_CONFIG_FILE = {PARM_BASE}/met_config/STATAnalysisConfig_wrapped
MODEL1 = GFS
MODEL1_OBTYPE = ADPUPA
# stat_analysis job info
STAT_ANALYSIS_JOB_NAME = aggregate_stat
# if using -dump_row, put in JOBS_ARGS "-dump_row [dump_row_file]"
# if using -out_stat, put in JOBS_ARGS "-out_stat [out_stat_file]"
# METplus will fill in filename
STAT_ANALYSIS_JOB_ARGS = -out_line_type CTS -out_thresh ==1 -out_stat [out_stat_file]
MODEL_LIST = {MODEL1}
FCST_LEAD_LIST = 24
LINE_TYPE_LIST = MPR
GROUP_LIST_ITEMS = MODEL_LIST
LOOP_LIST_ITEMS = FCST_LEAD_LIST
MODEL1_STAT_ANALYSIS_LOOKIN_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/mpr/IBL
STAT_ANALYSIS_OUTPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking
MODEL1_STAT_ANALYSIS_OUT_STAT_TEMPLATE = {model?fmt=%s}_ERA_IBLS_{lead?fmt=%H%M%S}L_CTS_CNT.stat
###
# StatAnalysis(sanal_blocks) Settings
# https://metplus.readthedocs.io/en/latest/Users_Guide/wrappers.html#statanalysis
###
# Stat Analysis for the Blocks
[sanal_blocks]
VALID_TIME_FMT = %Y%m%d
VALID_BEG = 20001201
VALID_END = 20170228
STAT_ANALYSIS_CONFIG_FILE = {PARM_BASE}/met_config/STATAnalysisConfig_wrapped
MODEL1 = GFS
MODEL1_OBTYPE = ADPUPA
STAT_ANALYSIS_JOB_NAME = aggregate_stat
STAT_ANALYSIS_JOB_ARGS = -out_line_type CTS -out_thresh ==1 -out_stat [out_stat_file]
MODEL_LIST = {MODEL1}
FCST_LEAD_LIST = 24
LINE_TYPE_LIST = MPR
GROUP_LIST_ITEMS = MODEL_LIST
LOOP_LIST_ITEMS = FCST_LEAD_LIST
MODEL1_STAT_ANALYSIS_LOOKIN_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/mpr/Blocks
STAT_ANALYSIS_OUTPUT_DIR = {OUTPUT_BASE}/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking
MODEL1_STAT_ANALYSIS_OUT_STAT_TEMPLATE = {model?fmt=%s}_ERA_Blocks_{lead?fmt=%H%M%S}L_CTS.stat
MET Configuration
METplus sets environment variables based on the values in the METplus configuration file. These variables are referenced in the MET configuration file. 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 not controlled by an environment variable, you can add additional environment variables to be set only within the METplus environment using the [user_env_vars] section of the METplus configuration files. See the ‘User Defined Config’ section on the ‘System Configuration’ page of the METplus User’s Guide for more information.
See the following files for more information about the environment variables set in this configuration file.
parm/use_cases/met_tool_wrapper/RegridDataPlane/RegridDataPlane.py parm/use_cases/met_tool_wrapper/PCPCombine/PCPCOmbine_derive.py parm/use_cases/met_tool_wrapper/PCPCombine/PCPCOmbine_subtract.py parm/use_cases/met_tool_wrapper/StatAnalysis/StatAnalysis.py
Python Scripts
This use case uses Python scripts to perform the blocking calculation
parm/use_cases/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking/Blocking_driver.py: This script calls the requested steps in the blocking analysis for a forecast, observation, or both. The possible steps are computing CBLs, plotting CBLs, computing IBLs, plotting IBLs, computing GIBLs, computing blocks, and plotting blocks.
metcalcpy/contributed/blocking_weather_regime/Blocking.py: This script runs the requested steps, containing the code for computing CBLs, computing IBLs, computing GIBLs, and computing blocks. See the METcalcpy Blocking Calculation Script for more information.
metcalcpy/contributed/blocking_weather_regime/Blocking_WeatherRegime_util.py: This script contains functions used by both the blocking anwd weather regime analysis, including the code for determining which steps the user wants to run, and finding and reading the input files in the format from the output pre-processing steps. See the METcalcpy Utility script for more information.
#!/usr/bin/env python3
import sys
import os
import numpy as np
import datetime
import netCDF4
import warnings
from metcalcpy.contributed.blocking_weather_regime.Blocking import BlockingCalculation
from metcalcpy.contributed.blocking_weather_regime.Blocking_WeatherRegime_util import parse_steps, write_mpr_file
from metplotpy.contributed.blocking_s2s import plot_blocking as pb
from metplotpy.contributed.blocking_s2s.CBL_plot import create_cbl_plot
def main():
steps_list_fcst,steps_list_obs = parse_steps()
if not steps_list_obs and not steps_list_fcst:
warnings.warn('No processing steps requested for either the model or observations,')
warnings.warn(' nothing will be run')
warnings.warn('Set FCST_STEPS and/or OBS_STEPS in the [user_env_vars] section to process data')
######################################################################
# Blocking Calculation and Plotting
######################################################################
# Set up the data
steps_fcst = BlockingCalculation('FCST')
steps_obs = BlockingCalculation('OBS')
# Check to see if there is a plot directory
oplot_dir = os.environ.get('BLOCKING_PLOT_OUTPUT_DIR','')
if not oplot_dir:
obase = os.environ['SCRIPT_OUTPUT_BASE']
oplot_dir = os.path.join(obase,'plots')
if not os.path.exists(oplot_dir):
os.makedirs(oplot_dir)
# Check to see if there is a mpr output directory
mpr_dir = os.environ.get('BLOCKING_MPR_OUTPUT_DIR','')
if not mpr_dir:
obase = os.environ['SCRIPT_OUTPUT_BASE']
mpr_dir = os.path.join(obase,'mpr')
# Check to see if CBL's are used from an obs climatology
use_cbl_obs = os.environ.get('USE_CBL_OBS','False').lower()
# Get the days per season
dseasons = int(os.environ['DAYS_PER_SEASON'])
# Grab the Anomaly (CBL) text files
obs_cbl_filetxt = os.environ.get('METPLUS_FILELIST_OBS_CBL_INPUT','')
fcst_cbl_filetxt = os.environ.get('METPLUS_FILELIST_FCST_CBL_INPUT','')
# Grab the Daily (IBL) text files
obs_ibl_filetxt = os.environ.get('METPLUS_FILELIST_OBS_IBL_INPUT','')
fcst_ibl_filetxt = os.environ.get('METPLUS_FILELIST_FCST_IBL_INPUT','')
# Calculate Central Blocking Latitude
if ("CBL" in steps_list_obs):
print('Computing Obs CBLs')
# Read in the list of CBL files
cbl_nseasons = int(os.environ['CBL_NUM_SEASONS'])
with open(obs_cbl_filetxt) as ocl:
obs_infiles = ocl.read().splitlines()
if (obs_infiles[0] == 'file_list'):
obs_infiles = obs_infiles[1:]
if len(obs_infiles) != (cbl_nseasons*dseasons):
raise Exception('Invalid Obs data; each year must contain the same date range to calculate seasonal averages.')
cbls_obs,lats_obs,lons_obs,mhweight_obs,cbl_time_obs = steps_obs.run_CBL(obs_infiles,cbl_nseasons,dseasons)
if ("CBL" in steps_list_fcst) and (use_cbl_obs == 'false'):
# Add in step to use obs for CBLS
print('Computing Forecast CBLs')
cbl_nseasons = int(os.environ['CBL_NUM_SEASONS'])
with open(fcst_cbl_filetxt) as fcl:
fcst_infiles = fcl.read().splitlines()
if (fcst_infiles[0] == 'file_list'):
fcst_infiles = fcst_infiles[1:]
if len(fcst_infiles) != (cbl_nseasons*dseasons):
raise Exception('Invalid Fcst data; each year must contain the same date range to calculate seasonal averages.')
cbls_fcst,lats_fcst,lons_fcst,mhweight_fcst,cbl_time_fcst = steps_fcst.run_CBL(fcst_infiles,cbl_nseasons,dseasons)
elif ("CBL" in steps_list_fcst) and (use_cbl_obs == 'true'):
if not ("CBL" in steps_list_obs):
raise Exception('Must run observed CBLs before using them as a forecast.')
cbls_fcst = cbls_obs
lats_fcst = lats_obs
lons_fcst = lons_obs
mhweight_fcst = mhweight_obs
cbl_time_fcst = cbl_time_obs
#Plot Central Blocking Latitude
if ("PLOTCBL" in steps_list_obs):
if not ("CBL" in steps_list_obs):
raise Exception('Must run observed CBLs before plotting them.')
print('Plotting Obs CBLs')
cbl_plot_mthstr = os.environ['OBS_CBL_PLOT_MTHSTR']
cbl_plot_outname = os.path.join(oplot_dir,os.environ.get('OBS_CBL_PLOT_OUTPUT_NAME','obs_cbl_avg'))
create_cbl_plot(lons_obs, lats_obs, cbls_obs, mhweight_obs, cbl_plot_mthstr, cbl_plot_outname,
do_averaging=True)
if ("PLOTCBL" in steps_list_fcst):
if not ("CBL" in steps_list_fcst):
raise Exception('Must run forecast CBLs before plotting them.')
print('Plotting Forecast CBLs')
cbl_plot_mthstr = os.environ['FCST_CBL_PLOT_MTHSTR']
cbl_plot_outname = os.path.join(oplot_dir,os.environ.get('FCST_CBL_PLOT_OUTPUT_NAME','fcst_cbl_avg'))
create_cbl_plot(lons_fcst, lats_fcst, cbls_fcst, mhweight_fcst, cbl_plot_mthstr, cbl_plot_outname,
do_averaging=True)
# Run IBL
if ("IBL" in steps_list_obs):
if not ("CBL" in steps_list_obs):
raise Exception('Must run observed CBLs before running IBLs.')
print('Computing Obs IBLs')
ibl_nseasons = int(os.environ['IBL_NUM_SEASONS'])
with open(obs_ibl_filetxt) as oil:
obs_infiles = oil.read().splitlines()
if (obs_infiles[0] == 'file_list'):
obs_infiles = obs_infiles[1:]
if len(obs_infiles) != (ibl_nseasons*dseasons):
raise Exception('Invalid Obs data; each year must contain the same date range to calculate seasonal averages.')
ibls_obs,ibl_time_obs = steps_obs.run_Calc_IBL(cbls_obs,obs_infiles,ibl_nseasons,dseasons)
daynum_obs = np.arange(0,len(ibls_obs[0,:,0]),1)
if ("IBL" in steps_list_fcst):
if (not "CBL" in steps_list_fcst):
raise Exception('Must run forecast CBLs or use observed CBLs before running IBLs.')
print('Computing Forecast IBLs')
ibl_nseasons = int(os.environ['IBL_NUM_SEASONS'])
with open(fcst_ibl_filetxt) as fil:
fcst_infiles = fil.read().splitlines()
if (fcst_infiles[0] == 'file_list'):
fcst_infiles = fcst_infiles[1:]
if len(fcst_infiles) != (ibl_nseasons*dseasons):
raise Exception('Invalid Fcst data; each year must contain the same date range to calculate seasonal averages.')
ibls_fcst,ibl_time_fcst = steps_fcst.run_Calc_IBL(cbls_fcst,fcst_infiles,ibl_nseasons,dseasons)
daynum_fcst = np.arange(0,len(ibls_fcst[0,:,0]),1)
if ("IBL" in steps_list_obs) and ("IBL" in steps_list_fcst):
# Print IBLs to output matched pair file
i_mpr_outdir = os.path.join(mpr_dir,'IBL')
if not os.path.exists(i_mpr_outdir):
os.makedirs(i_mpr_outdir)
modname = os.environ.get('MODEL_NAME','GFS')
maskname = os.environ.get('MASK_NAME','FULL')
ibl_outfile_prefix = os.path.join(i_mpr_outdir,'IBL_stat_'+modname)
cbls_avg = np.nanmean(cbls_obs,axis=0)
write_mpr_file(ibls_obs,ibls_fcst,cbls_avg,lons_obs,ibl_time_obs,ibl_time_fcst,modname,
'NA','IBLs','block','Z500','IBLs','block','Z500',maskname,'500',ibl_outfile_prefix)
# Plot IBLS
if("PLOTIBL" in steps_list_obs) and not ("PLOTIBL" in steps_list_fcst):
if not ("IBL" in steps_list_obs):
raise Exception('Must run observed IBLs before plotting them.')
print('Plotting Obs IBLs')
ibl_plot_title = os.environ.get('OBS_IBL_PLOT_TITLE','Instantaneous Blocked Longitude')
ibl_plot_outname = os.path.join(oplot_dir,os.environ.get('OBS_IBL_PLOT_OUTPUT_NAME','obs_IBL_Freq'))
ibl_plot_label1 = os.environ.get('IBL_PLOT_OBS_LABEL','')
pb.plot_ibls(ibls_obs,lons_obs,ibl_plot_title,ibl_plot_outname,label1=ibl_plot_label1)
elif ("PLOTIBL" in steps_list_fcst) and not ("PLOTIBL" in steps_list_obs):
if not ("IBL" in steps_list_fcst):
raise Exception('Must run forecast IBLs before plotting them.')
print('Plotting Forecast IBLs')
ibl_plot_title = os.environ.get('FCST_IBL_PLOT_TITLE','Instantaneous Blocked Longitude')
ibl_plot_outname = os.path.join(oplot_dir,os.environ.get('FCST_IBL_PLOT_OUTPUT_NAME','fcst_IBL_Freq'))
ibl_plot_label1 = os.environ.get('IBL_PLOT_FCST_LABEL','')
pb.plot_ibls(ibls_fcst,lons_fcst,ibl_plot_title,ibl_plot_outname,label1=ibl_plot_label1)
elif ("PLOTIBL" in steps_list_obs) and ("PLOTIBL" in steps_list_fcst):
if (not "IBL" in steps_list_obs) and (not "IBL" in steps_list_fcst):
raise Exception('Must run forecast and observed IBLs before plotting them.')
print('Plotting Obs and Forecast IBLs')
ibl_plot_title = os.environ['IBL_PLOT_TITLE']
ibl_plot_outname = os.path.join(oplot_dir,os.environ.get('IBL_PLOT_OUTPUT_NAME','IBL_Freq'))
#Check to see if there are plot legend labels
ibl_plot_label1 = os.environ.get('IBL_PLOT_OBS_LABEL','Observation')
ibl_plot_label2 = os.environ.get('IBL_PLOT_FCST_LABEL','Forecast')
pb.plot_ibls(ibls_obs,lons_obs,ibl_plot_title,ibl_plot_outname,data2=ibls_fcst,lon2=lons_fcst,
label1=ibl_plot_label1,label2=ibl_plot_label2)
# Run GIBL
if ("GIBL" in steps_list_obs):
if not ("IBL" in steps_list_obs):
raise Exception('Must run observed IBLs before running GIBLs.')
print('Computing Obs GIBLs')
gibls_obs = steps_obs.run_Calc_GIBL(ibls_obs,lons_obs)
if ("GIBL" in steps_list_fcst):
if not ("IBL" in steps_list_fcst):
raise Exception('Must run Forecast IBLs before running GIBLs.')
print('Computing Forecast GIBLs')
gibls_fcst = steps_fcst.run_Calc_GIBL(ibls_fcst,lons_fcst)
# Calc Blocks
if ("CALCBLOCKS" in steps_list_obs):
if not ("GIBL" in steps_list_obs):
raise Exception('Must run observed GIBLs before calculating blocks.')
print('Computing Obs Blocks')
block_freq_obs = steps_obs.run_Calc_Blocks(ibls_obs,gibls_obs,lons_obs,daynum_obs)
if ("CALCBLOCKS" in steps_list_fcst):
if not ("GIBL" in steps_list_fcst):
raise Exception('Must run Forecast GIBLs before calculating blocks.')
print('Computing Forecast Blocks')
block_freq_fcst = steps_fcst.run_Calc_Blocks(ibls_fcst,gibls_fcst,lons_fcst,daynum_fcst)
# Write out a Blocking MPR file if both obs and forecast blocking calculation performed
if ("CALCBLOCKS" in steps_list_obs) and ("CALCBLOCKS" in steps_list_fcst):
b_mpr_outdir = os.path.join(mpr_dir,'Blocks')
if not os.path.exists(b_mpr_outdir):
os.makedirs(b_mpr_outdir)
# Print Blocks to output matched pair file
modname = os.environ.get('MODEL_NAME','GFS')
maskname = os.environ.get('MASK_NAME','FULL')
blocks_outfile_prefix = os.path.join(b_mpr_outdir,'blocking_stat_'+modname)
cbls_avg = np.nanmean(cbls_obs,axis=0)
write_mpr_file(block_freq_obs,block_freq_fcst,cbls_avg,lons_obs,ibl_time_obs,ibl_time_fcst,modname,
'NA','Blocks','block','Z500','Blocks','block','Z500',maskname,'500',blocks_outfile_prefix)
# Plot Blocking Frequency
if ("PLOTBLOCKS" in steps_list_obs):
if not ("CALCBLOCKS" in steps_list_obs):
raise Exception('Must compute observed blocks before plotting them.')
print('Plotting Obs Blocks')
blocking_plot_title = os.environ.get('OBS_BLOCKING_PLOT_TITLE','Obs Blocking Frequency')
blocking_plot_outname = os.path.join(oplot_dir,os.environ.get('OBS_BLOCKING_PLOT_OUTPUT_NAME','obs_Block_Freq'))
pb.plot_blocks(block_freq_obs,gibls_obs,ibls_obs,lons_obs,blocking_plot_title,blocking_plot_outname)
if ("PLOTBLOCKS" in steps_list_fcst):
if not ("CALCBLOCKS" in steps_list_fcst):
raise Exception('Must compute forecast blocks before plotting them.')
print('Plotting Forecast Blocks')
blocking_plot_title = os.environ.get('FCST_BLOCKING_PLOT_TITLE','Forecast Blocking Frequency')
blocking_plot_outname = os.path.join(oplot_dir,os.environ.get('FCST_BLOCKING_PLOT_OUTPUT_NAME','fcst_Block_Freq'))
pb.plot_blocks(block_freq_fcst,gibls_fcst,ibls_fcst,lons_fcst,blocking_plot_title,blocking_plot_outname)
if __name__ == "__main__":
main()
Running METplus
This use case is run in the following ways:
Passing in UserScript_fcstGFS_obsERA_Blocking.py then a user-specific system configuration file:
run_metplus.py -c /path/to/METplus/parm/use_cases/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking.py -c /path/to/user_system.conf
Modifying the configurations in parm/metplus_config, then passing in UserScript_fcstGFS_obsERA_Blocking.py:
run_metplus.py -c /path/to/METplus/parm/use_cases/model_applications/s2s_mid_lat/UserScript_fcstGFS_obsERA_Blocking.py
The following variables must be set correctly:
INPUT_BASE - Path to directory where sample data tarballs are unpacked (See Datasets section to obtain tarballs). This is not required to run METplus, but it is required to run the examples in parm/use_cases
OUTPUT_BASE - Path where METplus output will be written. This must be in a location where you have write permissions
MET_INSTALL_DIR - Path to location where MET is installed locally
Example User Configuration File:
[dir]
INPUT_BASE = /path/to/sample/input/data
OUTPUT_BASE = /path/to/output/dir
MET_INSTALL_DIR = /path/to/met-X.Y
Expected Output
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 model_applications/s2s_mid_lat/Blocking (relative to OUTPUT_BASE) and will contain output for the steps requested. This may include the regridded data, daily averaged files, running mean files, and anomaly files. In addition, output CBL, IBL, and Blocking frequency plots can be generated. The location of these output plots can be specified as BLOCKING_PLOT_OUTPUT_DIR. If it is not specified, plots will be sent to OUTPUT_BASE/plots. MET format matched pair output will also be generated for IBLs and blocks if a user runs these steps on both the model and observation data. The location the matched pair output can be specified as BLOCKING_MPR_OUTPUT_DIR. If it is not specified, plots will be sent to OUTPUT_BASE/mpr. An output contingency table statistics line from stat_analysis is also generated from the IBL and Blocks matched pair files. The location of the output is set as STAT_ANALYSIS_OUTPUT_DIR.
Keywords
Note
RegridDataPlaneToolUseCase
PCPCombineToolUseCase
StatAnalysisToolUseCase
S2SAppUseCase
S2SMidLatAppUseCase
NetCDFFileUseCase
GRIB2FileUseCase
Navigate to the METplus Quick Search for Use Cases page to discover other similar use cases.
sphinx_gallery_thumbnail_path = ‘_static/s2s_mid_lat-UserScript_fcstGFS_obsERA_Blocking.png’
Total running time of the script: (0 minutes 0.000 seconds)