Feature #3024 and #3030 Series-Analysis GRAD

data/config/SeriesAnalysisConfig_default

@@ -125,6 +125,19 @@ mask = {
    poly = "";
 }
 
+////////////////////////////////////////////////////////////////////////////////
+
+//
+// Gradient statistics
+// May be set separately in each "obs.field" entry
+//
+gradient = {
+   dx = [ 1 ];
+   dy = [ 1 ];
+}
+
+////////////////////////////////////////////////////////////////////////////////
+
 //
 // Number of grid points to be processed concurrently. Set smaller to use less
 // memory but increase the number of passes through the data. If set <= 0, all
@@ -155,6 +168,7 @@ output_stats = {
    pstd   = [];
    pjc    = [];
    prc    = [];
+   grad   = [];
 }
 
 ////////////////////////////////////////////////////////////////////////////////

data/table_files/met_header_columns_V12.1.txt

@@ -11,7 +11,7 @@ V12.1 : STAT : SEEPS_MPR: VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID
 V12.1 : STAT : NBRCNT   : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL FBS FBS_BCL FBS_BCU FSS FSS_BCL FSS_BCU AFSS AFSS_BCL AFSS_BCU UFSS UFSS_BCL UFSS_BCU F_RATE F_RATE_BCL F_RATE_BCU O_RATE O_RATE_BCL O_RATE_BCU
 V12.1 : STAT : NBRCTC   : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL FY_OY FY_ON FN_OY FN_ON
 V12.1 : STAT : NBRCTS   : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL BASER BASER_NCL BASER_NCU BASER_BCL BASER_BCU FMEAN FMEAN_NCL FMEAN_NCU FMEAN_BCL FMEAN_BCU ACC ACC_NCL ACC_NCU ACC_BCL ACC_BCU FBIAS FBIAS_BCL FBIAS_BCU PODY PODY_NCL PODY_NCU PODY_BCL PODY_BCU PODN PODN_NCL PODN_NCU PODN_BCL PODN_BCU POFD POFD_NCL POFD_NCU POFD_BCL POFD_BCU FAR FAR_NCL FAR_NCU FAR_BCL FAR_BCU CSI CSI_NCL CSI_NCU CSI_BCL CSI_BCU GSS GSS_BCL GSS_BCU HK HK_NCL HK_NCU HK_BCL HK_BCU HSS HSS_BCL HSS_BCU ODDS ODDS_NCL ODDS_NCU ODDS_BCL ODDS_BCU LODDS LODDS_NCL LODDS_NCU LODDS_BCL LODDS_BCU ORSS ORSS_NCL ORSS_NCU ORSS_BCL ORSS_BCU EDS EDS_NCL EDS_NCU EDS_BCL EDS_BCU SEDS SEDS_NCL SEDS_NCU SEDS_BCL SEDS_BCU EDI EDI_NCL EDI_NCU EDI_BCL EDI_BCU SEDI SEDI_NCL SEDI_NCU SEDI_BCL SEDI_BCU BAGSS BAGSS_BCL BAGSS_BCU
-V12.1 : STAT : GRAD     : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL FGBAR OGBAR MGBAR EGBAR S1 S1_OG FGOG_RATIO DX DY
+V12.1 : STAT : GRAD     : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL FGBAR OGBAR MGBAR EGBAR S1 S1_OG FGOG_RATIO DX DY FGMAG OGMAG MAG_RMSE LAPLACE_RMSE
 V12.1 : STAT : DMAP     : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL FY OY FBIAS BADDELEY HAUSDORFF MED_FO MED_OF MED_MIN MED_MAX MED_MEAN FOM_FO FOM_OF FOM_MIN FOM_MAX FOM_MEAN ZHU_FO ZHU_OF ZHU_MIN ZHU_MAX ZHU_MEAN G GBETA BETA_VALUE 
 V12.1 : STAT : ORANK    : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL INDEX OBS_SID OBS_LAT OBS_LON OBS_LVL OBS_ELV OBS PIT RANK N_ENS_VLD (N_ENS) ENS_[0-9]* OBS_QC ENS_MEAN OBS_CLIMO_MEAN SPREAD ENS_MEAN_OERR SPREAD_OERR SPREAD_PLUS_OERR OBS_CLIMO_STDEV FCST_CLIMO_MEAN FCST_CLIMO_STDEV
 V12.1 : STAT : PCT      : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL (N_THRESH) THRESH_[0-9]* OY_[0-9]* ON_[0-9]*

data/table_files/stat_column_description.txt

@@ -272,3 +272,17 @@ PRC_N_THRESH             "Number of probability thresholds"
 PRC_THRESH_I             "The i-th probability threshold"
 PRC_PODY_I               "Probability of detecting yes when forecast is between the i-th and (i+1)-th probability thresholds"
 PRC_POFD_I               "Probability of false detection when forecast is between the i-th and (i+1)-th probability thresholds"
+GRAD_TOTAL               "Total number of matched pairs"
+GRAD_FGBAR               "Mean of absolute value of forecast gradients"
+GRAD_OGBAR               "Mean of absolute value of observed gradients"
+GRAD_MGBAR               "Mean of maximum of absolute values of forecast and observed gradients"
+GRAD_EGBAR               "Mean of absolute value of forecast minus observed gradients"
+GRAD_S1                  "S1 score"
+GRAD_S1_OG               "S1 score with respect to observed gradient"
+GRAD_FGOG_RATIO          "Ratio of forecast and observed gradients"
+GRAD_DX                  "Gradient size in the X-direction"
+GRAD_DY                  "Gradient size in the Y-direction"
+GRAD_FGMAG               "Magnitude of the forecast gradient when the X and Y-directions are interpreted as a vector"
+GRAD_OGMAG               "Magnitude of the observed gradient when the X and Y-directions are interpreted as a vector"
+GRAD_MAG_RMSE            "Root mean squared difference of the forecast gradient magnitude minus the observed gradient magnitude"
+GRAD_LAPLACE_RMSE        "Root mean squared difference of the sum of the forecast X and Y-gradients minus the number of the observed X and Y-gradients"

docs/Users_Guide/appendixC.rst

@@ -739,7 +739,7 @@ These statistics require climatological values for the wind vector components, w
 Gradient Values
 ---------------
 
-Called "TOTAL", "FGBAR", "OGBAR", "MGBAR", "EGBAR", "S1", "S1_OG", and "FGOG_RATIO" in GRAD output :numref:`table_GS_format_info_GRAD`
+Called "TOTAL", "FGBAR", "OGBAR", "MGBAR", "EGBAR", "S1", "S1_OG", "FGOG_RATIO", "FGMAG", "OGMAG", "MAG_RMSE", and "LAPLACE_RMSE" in GRAD output :numref:`table_GS_format_info_GRAD`
 
 These statistics are only computed by the Grid-Stat tool and require vectors. Here :math:`\nabla` is the gradient operator, which in this applications signifies the difference between adjacent grid points in both the grid-x and grid-y directions. TOTAL is the count of grid locations used in the calculations. The remaining measures are defined below:
 
@@ -778,6 +778,39 @@ where the weights are applied at each grid location, with values assigned accord
 
     \text{FGOG_RATIO} = \frac{\text{FGBAR}}{\text{OGBAR}}
 
+The following statistics are computed using the magnitude of the vectors formed by the forecast and observation gradients in the grid-x and grid-y directions.
+
+.. math::
+  \text{FGMAG} = \text{Mean}(|| \nabla f ||) = \frac{1}{n} \sum_{i=1}^n \sqrt{\nabla {f_x}_i^2 + \nabla {f_y}_i^2}
+
+  \text{OGMAG} = \text{Mean}(|| \nabla o ||) = \frac{1}{n} \sum_{i=1}^n \sqrt{\nabla {o_x}_i^2 + \nabla {o_y}_i^2}
+
+.. only:: latex
+
+  .. math::
+
+    \text{MAG\_RMSE} = \sqrt{ \frac{1}{n} \sum_{i=1}^n {(|| \nabla f_i || - || \nabla o_i ||)}^2 }
+
+.. only:: html
+
+  .. math::
+
+    \text{MAG_RMSE} = \sqrt{ \frac{1}{n} \sum_{i=1}^n {(|| \nabla f_i || - || \nabla o_i ||)}^2 }
+
+Laplace RMSE is very similar to gradient RMSE, but instead of taking the magnitude of the gradient vector at each
+point, we compute the divergence of the gradient.
+
+.. only:: latex
+
+  .. math::
+
+    \text{LAPLACE\_RMSE} = \sqrt{ \frac{1}{n} \sum_{i=1}^n { ((\nabla {f_x}_i + \nabla {f_y}_i) - (\nabla {o_x}_i + \nabla {o_y}_i))^2 }}
+
+.. only:: html
+
+  .. math::
+
+    \text{LAPLACE_RMSE} = \sqrt{ \frac{1}{n} \sum_{i=1}^n { ((\nabla {f_x}_i + \nabla {f_y}_i) - (\nabla {o_x}_i + \nabla {o_y}_i))^2 }}
 
 MET Verification Measures for Probabilistic Forecasts
 =====================================================

docs/Users_Guide/config_options.rst

@@ -2232,6 +2232,32 @@ This dictionary may include the following entries:
      prob_cat_thresh = [];
   }
 
+
+.. _gradient:
+
+gradient
+--------
+
+:ref:`gradient <gradient>`
+
+The "gradient" entry is a dictionary which specifies the number and size
+of gradients to be computed and applies to both Grid-Stat and
+Series-Analysis. The "dx" and "dy" entries specify the size of the
+gradients in grid units in the X and Y dimensions, respectively. "dx"
+and "dy" are arrays of integers (positive or negative) which must have the
+same length, and the GRAD output line type will be computed separately for
+each entry. When computing gradients, the value at the (x, y) grid point
+is replaced by the value at the (x+dx, y+dy) grid point minus the value at
+(x, y). This configuration option may be set separately in each "obs.field"
+entry.
+
+.. code-block:: none
+
+  gradient = {
+     dx = [ 1 ];
+     dy = [ 1 ];
+   }
+
 output_flag
 -----------
 

docs/Users_Guide/grid-stat.rst

@@ -98,6 +98,13 @@ The S1 score has been in historical use for verification of forecasts, particula
 
 Differences are computed in both of the horizontal grid directions and is not a true mathematical gradient. Because the S1 score focuses on differences only, any bias in the forecast will not be measured. Further, the score depends on the domain and spacing of the grid, so can only be compared on forecasts with identical grids.
 
+As described in :ref:`Ebert-Uphoff et al., 2024 <Ebert-Uphoff-2024>`, statistics based
+on the magnitude of the forecast and observed gradients are also provided. Similiar to
+the S1 score, the root-mean-squared error of the magnitude of the gradients and their
+divergence quantify the similarity in the texture of the fields, with 0 being a perfect
+score. These gradient-based statistics assess the difference in smoothness between the
+two fields but not the accuracy of the forecast.
+
 Distance Maps
 -------------
 
@@ -263,6 +270,7 @@ __________________________
   eclv_points      = 0.05;
   hss_ec_value     = NA;
   rank_corr_flag   = TRUE;
+  gradient         = { dx = [ 1 ]; dy = [ 1 ]; }
   grid_weight_flag = NONE;
   tmp_dir          = "/tmp";
   output_prefix    = "";
@@ -319,22 +327,6 @@ The available wave numbers start at 0 (the mean across each row of data) and end
 
 _____________________
 
-.. _gradient:
-
-:ref:`gradient <gradient>`
-
-.. code-block:: none
-
-  gradient = {
-     dx = [ 1 ];
-     dy = [ 1 ];
-   }
-
-
-The **gradient** entry is a dictionary which specifies the number and size of gradients to be computed. The **dx** and **dy** entries specify the size of the gradients in grid units in the X and Y dimensions, respectively. **dx** and **dy** are arrays of integers (positive or negative) which must have the same length, and the GRAD output line type will be computed separately for each entry. When computing gradients, the value at the (x, y) grid point is replaced by the value at the (x+dx, y+dy) grid point minus the value at (x, y). This configuration option may be set separately in each **obs.field** entry.
-
-____________________
-
 .. _distance_map:
 
 :ref:`distance_map <distance_map>`
@@ -832,11 +824,27 @@ The format of the STAT and ASCII output of the Grid-Stat tool are the same as th
   * - 33
     - DX
     - Gradient size in the X-direction
-    - Integer
+    - Integer 
   * - 34
     - DY
     - Gradient size in the Y-direction
-    - Integer
+    - Integer 
+  * - 35
+    - FGMAG 
+    - Magnitude of the forecast gradient when the X and Y-directions are interpreted as a vector
+    - Double
+  * - 36
+    - OGMAG 
+    - Magnitude of the observed gradient when the X and Y-directions are intrepreted as a vector
+    - Double
+  * - 37
+    - MAG_RMSE 
+    - Root mean squared difference of the forecast gradient magnitude minus the observed gradient magnitude
+    - Double
+  * - 38
+    - LAPLACE_RMSE 
+    - Root mean squared difference of the sum of the forecast X and Y-gradients minus the sum of the observed X and Y-gradients
+    - Double
 
 .. _table_GS_format_info_DMAP:
 

docs/Users_Guide/refs.rst

@@ -130,6 +130,12 @@ References
 |   a review and proposed framework. *Meteorological Applications*, 15, 51-64.
 |
 
+.. _Ebert-Uphoff-2024:
+
+| Ebert-Uphoff, I.,, 2024: An Investigation of Metrics to Evaluate the Sharpness
+|   in AI-Generated Meteorological Imagery. *Draft version - Jan 26, 2024*
+|
+
 .. _Eckel-2012:
 
 | Eckel, F. A., M.S. Allen, M. C. Sittel, 2012: Estimation of Ambiguity in
@@ -278,7 +284,7 @@ References
 
 .. _North-2022:
 
-| North, R.C.,  M.P. Mittermaier, S.F. Milton, 2022. *Using SEEPS with a*
+| North, R.C., M.P. Mittermaier, S.F. Milton, 2022. *Using SEEPS with a*
 |   TRMM-derived Climatology to Assess Global NWP Precipitation Forecast Skill.
 |   *Monthly Weather Review*, 150, 135-155.
 |   https://doi.org/10.1175/MWR-D-20-0347.1

docs/Users_Guide/series-analysis.rst

@@ -116,6 +116,7 @@ ____________________
   boot           = { interval = PCTILE; rep_prop = 1.0; n_rep = 1000;
                      rng = "mt19937"; seed = ""; }
   mask           = { grid = [ "FULL" ]; poly = []; }
+  gradient       = { dx = [ 1 ]; dy = [ 1 ]; }
   hss_ec_value   = NA;
   rank_corr_flag = TRUE;
   tmp_dir        = "/tmp";
@@ -139,7 +140,6 @@ ____________________
 
 Ratio of valid matched pairs for the series of values at each grid point required to compute statistics. Set to a lower proportion to allow some missing values. Setting it to 1.0 requires that every data point be valid over the series to compute statistics.
 
-
 ____________________
 
 .. code-block:: none
@@ -157,6 +157,7 @@ ____________________
      pstd   = [];
      pjc    = [];
      prc    = [];
+     grad   = [];
   }
 
 The output_stats array controls the type of output that the Series-Analysis tool generates. Each flag corresponds to an output line type in the STAT file and is used to specify the comma-separated list of statistics to be computed. Use the column names from the tables listed below to specify the statistics. The output flags correspond to the following types of output line types:
@@ -185,4 +186,6 @@ The output_stats array controls the type of output that the Series-Analysis tool
 
 12. PRC for Receiver Operating Characteristic for Probabilistic forecasts (See :numref:`table_PS_format_info_PRC`)
 
-.. note:: When the -input option is used, all partial sum and contingency table count columns are required to aggregate statistics across multiple runs. To facilitate this, the output_stats entries for the CTC, SL1L2, SAL1L2, and PCT line types can be set to "ALL" to indicate that all available columns for those line types should be written.
+13. GRAD for Gradient Statistics (See :numref:`table_GS_format_info_GRAD`)
+
+.. note:: When the -input option is used, all partial sum and contingency table count columns are required to aggregate statistics across multiple runs. To facilitate this, the output_stats entries for the CTC, SL1L2, SAL1L2, PCT, and GRAD line types can be set to "ALL" to indicate that all available columns for those line types should be written.

internal/test_unit/config/SeriesAnalysisConfig

@@ -109,6 +109,19 @@ mask = {
    poly = "${MASK_POLY}";
 }
 
+////////////////////////////////////////////////////////////////////////////////
+
+//
+// Gradient statistics
+// May be set separately in each "obs.field" entry
+//
+gradient = {
+   dx = [ 1, 3 ];
+   dy = [ 1, 3 ];
+}
+
+////////////////////////////////////////////////////////////////////////////////
+
 //
 // Number of grid points to be processed concurrently. Set smaller to use less
 // memory but increase the number of passes through the data. If set <= 0, all
@@ -139,6 +152,7 @@ output_stats = {
    pstd   = [ ${PSTD_STATS}   ];
    pjc    = [ ${PJC_STATS}    ];
    prc    = [ ${PRC_STATS}    ];
+   grad   = [ ${GRAD_STATS}   ];
 }
 
 ////////////////////////////////////////////////////////////////////////////////

internal/test_unit/config/SeriesAnalysisConfig_climo

@@ -113,6 +113,19 @@ mask = {
    poly = "";
 }
 
+////////////////////////////////////////////////////////////////////////////////
+
+//
+// Gradient statistics
+// May be set separately in each "obs.field" entry
+//
+gradient = {
+   dx = [ 1 ];
+   dy = [ 1 ];
+}
+
+////////////////////////////////////////////////////////////////////////////////
+
 //
 // Number of grid points to be processed concurrently. Set smaller to use less
 // memory but increase the number of passes through the data. If set <= 0, all
@@ -133,16 +146,17 @@ vld_thresh = 0.5;
 output_stats = {
    fho    = [ "TOTAL", "F_RATE", "H_RATE", "O_RATE" ];
    ctc    = [ "ALL" ];
-   cts    = [ ];
-   mctc   = [ ];
-   mcts   = [ ];
+   cts    = [];
+   mctc   = [];
+   mcts   = [];
    cnt    = [ "TOTAL", "RMSE", "ANOM_CORR", "RMSFA", "RMSOA" ];
    sl1l2  = [ "ALL" ];
    sal1l2 = [ "ALL" ];
-   pct    = [ ];
-   pstd   = [ ];
-   pjc    = [ ];
-   prc    = [ ];
+   pct    = [];
+   pstd   = [];
+   pjc    = [];
+   prc    = [];
+   grad   = [];
 }
 
 ////////////////////////////////////////////////////////////////////////////////

internal/test_unit/config/SeriesAnalysisConfig_climo_prob

@@ -122,6 +122,19 @@ mask = {
    poly = "";
 }
 
+////////////////////////////////////////////////////////////////////////////////
+
+//
+// Gradient statistics
+// May be set separately in each "obs.field" entry
+//
+gradient = {
+   dx = [ 1 ];
+   dy = [ 1 ];
+}
+
+////////////////////////////////////////////////////////////////////////////////
+
 //
 // Number of grid points to be processed concurrently. Set smaller to use less
 // memory but increase the number of passes through the data. If set <= 0, all
@@ -140,18 +153,19 @@ vld_thresh = 0.5;
 // Statistical output types
 //
 output_stats = {
-   fho    = [ ];
-   ctc    = [ ];
-   cts    = [ ];
-   mctc   = [ ];
-   mcts   = [ ];
-   cnt    = [ ];
-   sl1l2  = [ ];
-   sal1l2 = [ ];
+   fho    = [];
+   ctc    = [];
+   cts    = [];
+   mctc   = [];
+   mcts   = [];
+   cnt    = [];
+   sl1l2  = [];
+   sal1l2 = [];
    pct    = [ "ALL" ];
    pstd   = [ "TOTAL", "ROC_AUC", "BRIER", "BRIERCL", "BSS", "BSS_SMPL" ];
-   pjc    = [ ];
-   prc    = [ ];
+   pjc    = [];
+   prc    = [];
+   grad   = [];
 }
 
 ////////////////////////////////////////////////////////////////////////////////