Feature 2055 nbm grib2 update (#2234)

* Per issue #2055, in read_grib2_record_list() added new PDS numbers 6, 10, 15 (for NBM data). SL

* Per issue #2055: in read_grib2_record_list(), in section that validates PDS number changed, for undefined PDS numbers, changed Error to Warning and continue (instead of exit). SL

* Per issue #2055, created new entries for NBM grib2 data fields. SL

* Per issue #2055. Added more NBM fields based on a 12z, f003 file. SL

* Per issue #2055, after interogating several NBM lead-time files, added TMAX field. SL

* Per issue #2055, added grib2_nbm.txt to the list. SL

* Per #2055, update the grib2_nbm.txt file entries in the table_files Makefiles.

* Per issue #2055, fixed syntax entry for TSTM field. SL

* Per issue #2055, fixed a few more syntax errors for some of the fields. SL

* Per issue #2055, in read_grib2_record_list(), added section to set level values for template number 6 (NBM data). Also added some temporary print statements. SL

* Per issue  #2055, in read_grib2_record_list(), added code to set levels for table number 10 (like 6). SL

* Per #2055, adding support for new GRIB2_perc_val configuration option to filter records based on the forecast percentile value used by GRIB2 product definition templates 6 and 10. Seth, note that I also included the probability filtering logic tweak we discussed, since it was already present in that same file.

* Per issue #2055: in find_record_matches(), in section that resolves Prob fields, changed 'break' to 'continue' for non Probs. SL

* Per issue #2055, in read_grib2_record_list() cleaned up some commented out code. SL

* Per issue #2055, in read_grib2_record_list() cleaned up all print statements. SL

* Per issue #2055, added 5 new unit tests for NBM data. SL

* Update data2d_grib2.cc

* Update data2d_grib2.cc

Co-authored-by: Seth Linden <linden@seneca.rap.ucar.edu>
Co-authored-by: John Halley Gotway <johnhg@ucar.edu>

data/table_files/Makefile.am

@@ -92,7 +92,8 @@ tablefiles_DATA = \
 	grib2_bom.txt \
 	grib2_gpp_af.txt \
 	grib2_mrms.txt \
-	grib2_ndfd.txt
+	grib2_ndfd.txt \
+	grib2_nbm.txt
 
 EXTRA_DIST = ${tablefiles_DATA}
 

data/table_files/Makefile.in

@@ -368,7 +368,8 @@ tablefiles_DATA = \
 	grib2_bom.txt \
 	grib2_gpp_af.txt \
 	grib2_mrms.txt \
-	grib2_ndfd.txt
+	grib2_ndfd.txt \
+	grib2_nbm.txt
 
 EXTRA_DIST = ${tablefiles_DATA}
 MAINTAINERCLEANFILES = Makefile.in

data/table_files/grib2_nbm.txt

@@ -0,0 +1,45 @@
+GRIB2
+0  1 0 255 7 1 0  21   "APTMP"  "Apparent Temperature"  "K"
+0  1 0 255 7 1 7  6    "CAPE"  "Convective Available Potential Energy"  "J/kg"
+0  1 0 255 7 1 19 239  "CWASP" "Craven-Wiedenfeld Aggregate Severe Parameter" "Numeric"
+0  1 0 255 7 1 0  6    "DPT" "Dew Point Temperature" "K"
+0  1 0 255 7 1 4  7    "DSWRF"  "Downward Short-Wave Radiation Flux"  "W/m^2"
+0  1 0 255 7 1 19 238  "ELLINX"  "Ellrod Index" "non-dim"
+0  1 0 255 7 1 1  228  "FICEAC"  "Flat Ice Accumulation (FRAM)"  "kg/m^2"
+10 1 0 255 7 1 3  204  "FRZSPR"  "Freezing SprayFreezing Spray"  "non-dim"
+0  1 0 255 7 1 19 10  "TURB"  "Turbulence"  "non-dim"
+2  1 0 255 7 1 4  4  "FOSINDX"  "Fosberg Index" "Numeric"
+2  1 0 255 7 1 4  2  "HINDEX"  "Haines Index"  "Numeric"
+0  1 0 255 7 1 0  5  "TMIN"  "Minimum Temperature"  "K"
+0  1 0 255 7 1 0  4  "TMAX"  "Maximum Temperature"  "K"
+0  1 0 255 7 1 19 3  "MIXHT"  "Mixed Layer Depth"  "m"
+0  1 0 255 7 1 0  27  "WETBT"  "Wet Bulb Temperature"  "K"
+0  1 0 255 7 1 7  201  "BNEGELAY" "Bourgoiun Negative Energy Layer"  "J/kg"
+0  1 0 255 7 1 1  234  "PCPDUR"  "Precipitation Duration"  "hour"
+0  1 0 255 7 1 1  8  "APCP"  "Total Precipitation"  "kg/m^2"
+0  1 0 255 7 1 7  202 "BPOSELAY" "Bourgoiun Positive Energy Layer (2k ft AGL to 400 hPa)" "J/kg"
+0  1 0 255 7 1 1  53  "TSRWE"  "Total Snowfall Rate Water Equivalent" "prob"
+0  1 0 255 7 1 1  226 "PWTHER"  "Predominant Weather"  "Numeric"
+0  1 0 255 7 1 1  19  "PTYPE"  "Precipitation Type" "prob"
+0  1 0 255 7 1 1  232  "PROBCIP"  "Probability Cloud Ice Present" "prob"
+0  1 0 255 7 1 19 2  "TSTM"  "Thunderstorm Probability"  "%"
+0  1 0 255 7 1 1  1  "RH"  "Relative Humidity"  "%"
+10 1 0 255 7 1 0  3  "HTSGW"  "Significant Height of Combined Wind Waves and Swell"  "m"
+0  1 0 255 7 1 6  1  "TCDC"  "Total Cloud Cover"  "%"
+0  1 0 255 7 1 1  233  "SNOWLR"  "Snow Liquid ratio"  "kg/kg"
+0  1 0 255 7 1 1  29  "ASNOW"  "Total Snowfall"  "kg/m^2"
+0  1 0 255 7 1 19 236  "SNOWLVL"  "Snow Level"  "m"
+0  1 0 255 7 1 0  0  "TMP"  "Temperature"  "K"
+0  1 0 255 7 1 2  226  "TRWDIR"  "Transport Wind Direction"   "deg"
+0  1 0 255 7 1 2  225  "TRWSPD"  "Transport Wind Speed"  "m/s"
+0  1 0 255 7 1 2  224  "VRATE"  "Ventilation Rate"  "m^2/s"
+0  1 0 255 7 1 2  0  "WDIR"  "Wind Direction (from which blowing)"  "deg"
+0  1 0 255 7 1 2  22  "GUST"  "Wind Speed (Gust)"  "m/s"
+0  1 0 255 7 1 2  1  "WIND"  "Wind Speed"  "m/s"
+0  1 0 255 7 1 6  13  "CEIL"  "Ceiling"  "m"
+0  1 0 255 7 1 19 237  "DRYTPROB"  "Dry Thunderstorm Probability"  "prob"
+0  1 0 255 7 1 16 3  "RETOP"  "Echo Top"  "m"
+0  1 0 255 7 1 3  6  "DIST"  "Geometric Height"  "m"
+0  1 0 255 7 1 16 198  "MAXREF"  "Hourly Maximum of Simulated Reflectivity at 1 km AGL"  "dB"
+0  1 0 255 7 1 15 3  "VIL"  "Vertically-Integrated Liquid Water"  "kg/m^2"
+0  1 0 255 7 1 19 0  "VIS"  "Visibility"  "m"

docs/Users_Guide/config_options.rst

@@ -974,6 +974,10 @@ File-format specific settings for the "field" entry:
 
     * The "GRIB2_stat_type" is an integer specifying the statistical
       processing type (Table 4.10).
+
+    * The "GRIB2_perc_val" is an integer specifying the requested percentile
+      value (0 to 100) to be used. This applies only to GRIB2 product
+      definition templates 4.6 and 4.10.
 
     * The "GRIB2_ipdtmpl_index" and "GRIB2_ipdtmpl_val" entries are arrays
       of integers which specify the product description template values to

internal/test_unit/xml/unit_plot_data_plane.xml

@@ -437,5 +437,75 @@
       <ps>&OUTPUT_DIR;/plot_data_plane/CESM_b.e15.B1850.f09_g16.pi_control.all_ga7cpl24.66.cam.i.0083-01-01-00000.ps</ps>
     </output>
   </test>
+
+  <test name="plot_data_plane_GRIB2_NBM_CWASP_L0">
+    <exec>&MET_BIN;/plot_data_plane</exec>
+    <param> \
+      &DATA_DIR_MODEL;/grib2/nbm/blend.t13z.core.f119.co.grib2 \
+      &OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_CWASP_L0.ps \
+      'name="CWASP"; level="L0"; GRIB2_pdt=2;' \
+      -title "GRIB2 NBM Craven-Wiedenfeld Aggregate Severe Parameter at surface" \
+      -v 1
+    </param>
+    <output>
+      <ps>&OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_CWASP_L0.ps</ps>
+    </output>
+  </test>
+
+  <test name="plot_data_plane_GRIB2_NBM_CWASP_PERC_5">
+    <exec>&MET_BIN;/plot_data_plane</exec>
+    <param> \
+      &DATA_DIR_MODEL;/grib2/nbm/blend.t13z.core.f119.co.grib2 \
+      &OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_CWASP_Perc-5.ps \
+      'name="CWASP"; level="L0"; GRIB2_perc_val=5;' \
+      -title "GRIB2 NBM Craven-Wiedenfeld Aggregate Severe Parameter 5% Level" \
+      -v 1
+    </param>
+    <output>
+      <ps>&OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_CWASP_Perc-5.ps</ps>
+    </output>
+  </test>
+
+  <test name="plot_data_plane_GRIB2_NBM_CWASP_PROB_50">
+    <exec>&MET_BIN;/plot_data_plane</exec>
+    <param> \
+      &DATA_DIR_MODEL;/grib2/nbm/blend.t13z.core.f119.co.grib2 \
+      &OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_CWASP_Prob-50.ps \
+      'name = "PROB"; level="L0"; prob = { name = "CWASP"; thresh_lo = 50; }' \
+      -title "GRIB2 NBM Craven-Wiedenfeld Aggregate Severe Parameter Prob > 50" \
+      -v 1
+    </param>
+    <output>
+      <ps>&OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_CWASP_Prob-50.ps</ps>
+    </output>
+  </test>
 
+  <test name="plot_data_plane_GRIB2_NBM_WETBT_MIXED_LEVELS">
+    <exec>&MET_BIN;/plot_data_plane</exec>
+    <param> \
+      &DATA_DIR_MODEL;/grib2/nbm/blend.t13z.core.f119.co.grib2 \
+      &OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_WETBT_L610-40000.ps \
+      'name="WETBT"; level="L610-40000";' \
+      -title "GRIB2 NBM Wet Bulb Temperature 610m above ground - 400mb" \
+      -v 1
+    </param>
+    <output>
+      <ps>&OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_WETBT_L610-40000.ps</ps>
+    </output>
+  </test>
+
+  <test name="plot_data_plane_GRIB2_NBM_FICEAC_A48_PERC_10">
+    <exec>&MET_BIN;/plot_data_plane</exec>
+    <param> \
+      &DATA_DIR_MODEL;/grib2/nbm/blend.t13z.core.f119.co.grib2 \
+      &OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_FICEAC_A48_Perc-10.ps \
+      'name="FICEAC"; level="A48"; GRIB2_perc_val=10;' \
+      -title "GRIB2 NBM 48-hour Flat Ice Accumulation (FRAM) 10% Level" \
+      -v 1
+    </param>
+    <output>
+      <ps>&OUTPUT_DIR;/plot_data_plane/nbm_2022021513_F119_GRIB2_FICEAC_A48_Perc-10.ps</ps>
+    </output>
+  </test>
+
 </met_test>

src/basic/vx_config/config_constants.h

@@ -525,6 +525,7 @@ static const char conf_key_GRIB2_process[]     = "GRIB2_process";
 static const char conf_key_GRIB2_ens_type[]    = "GRIB2_ens_type";
 static const char conf_key_GRIB2_der_type[]    = "GRIB2_der_type";
 static const char conf_key_GRIB2_stat_type[]   = "GRIB2_stat_type";
+static const char conf_key_GRIB2_perc_val[]    = "GRIB2_perc_val";
 static const char conf_key_GRIB2_ipdtmpl_index[] = "GRIB2_ipdtmpl_index";
 static const char conf_key_GRIB2_ipdtmpl_val[]   = "GRIB2_ipdtmpl_val";
 static const char conf_key_level[]             = "level";

src/libcode/vx_data2d_grib2/data2d_grib2.cc

@@ -432,7 +432,8 @@ void MetGrib2DataFile::find_record_matches( VarInfoGrib2* vinfo,
           (!is_bad_data(vinfo->process())   && vinfo->process()   != (*it)->Process  ) ||
           (!is_bad_data(vinfo->ens_type())  && vinfo->ens_type()  != (*it)->EnsType  ) ||
           (!is_bad_data(vinfo->der_type())  && vinfo->der_type()  != (*it)->DerType  ) ||
-          (!is_bad_data(vinfo->stat_type()) && vinfo->stat_type() != (*it)->StatType ) ){
+          (!is_bad_data(vinfo->stat_type()) && vinfo->stat_type() != (*it)->StatType ) ||
+          (!is_bad_data(vinfo->perc_val())  && vinfo->perc_val()  != (*it)->PercVal  ) ){
          continue;
       }
 
@@ -523,10 +524,13 @@ void MetGrib2DataFile::find_record_matches( VarInfoGrib2* vinfo,
          }  //  END: if( level match )
 
          //  if seeking a probabilistic field, check the prob info
-         if( (rec_match_ex || rec_match_rn) && vinfo->p_flag() && (*it)->ProbFlag ){
+         if( (rec_match_ex || rec_match_rn) && vinfo->p_flag() ) {
 
             rec_match_ex = rec_match_rn = false;
-
+
+            //  no match unless the data contains probabilities
+            if( !(*it)->ProbFlag ) { continue; }
+
             SingleThresh v_thr_lo = vinfo->p_thresh_lo();
             SingleThresh v_thr_hi = vinfo->p_thresh_hi();
 
@@ -695,17 +699,24 @@ void MetGrib2DataFile::read_grib2_record_list() {
               1 != gfld->ipdtnum &&     //  individual ensemble forecast, control and perturbed, at a horizontal level or in a horizontal layer at a point in time
               2 != gfld->ipdtnum &&     //  ensemble mean
               5 != gfld->ipdtnum &&     //  probability forecast
+              6 != gfld->ipdtnum &&     //  percentile forecasts at a horizontal level or in a horizontal layer at a point in time
               8 != gfld->ipdtnum &&     //  accumulation forecast
               9 != gfld->ipdtnum &&     //  probabilistic accumulation forecast
+             10 != gfld->ipdtnum &&     //  percentile forecasts at a horizontal level or in a horizontal layer in a continuous or non-continuous time interval
              11 != gfld->ipdtnum &&     //  individual ensemble forecast, control and perturbed, at a horizontal level or in a horizontal layer, in a continuous or non-continuous time interval
              12 != gfld->ipdtnum &&     //  derived accumulation forecast (?)
+             15 != gfld->ipdtnum &&     //  Average, accumulation, extreme values or other statistically-processed values over a spatial area at a horizontal level or in a horizontal layer at a point in time
              46 != gfld->ipdtnum &&     //  average, accumulation, and/or extreme values or other statistically processed values at a horizontal level or in a horizontal layer in a continuous or non-continuous time interval for aerosol.
              48 != gfld->ipdtnum ){     //  analysis or forecast at a horizontal level or in a horizontal layer at a point in time for aerosol.
-            mlog << Error << "\nMetGrib2DataFile::data_plane() -> "
+
+            // Print Warning, continue
+            mlog << Warning << "\nMetGrib2DataFile::data_plane() -> "
                  << "PDS template number ("
                  << gfld->ipdtnum << ") is not supported. "
-                 << "Please create a new post with this information in the METplus GitHub Discussions forum at https://github.com/dtcenter/METplus/discussions\n\n";
-            exit(1);
+                 << "Please create a new post with this information in the METplus GitHub Discussions forum at https://github.com/dtcenter/METplus/discussions. "
+                 << "Continuing to try and get the record information.\n\n";
+
+            continue;
          }
 
          //  store the record information
@@ -729,15 +740,16 @@ void MetGrib2DataFile::read_grib2_record_list() {
          }
 
          //  store the full pdtmpl values
-         for(int j=0; j < gfld->ipdtlen; j++){ rec->IPDTmpl.add((int) gfld->ipdtmpl[j]); }
-
+         for(int j=0; j < gfld->ipdtlen; j++) {
+            rec->IPDTmpl.add((int) gfld->ipdtmpl[j]);
+         }
+
          //  check for template number 46
          if( gfld->ipdtnum == 46 ) {
             rec->LvlVal1 = scaled2dbl(gfld->ipdtmpl[16], gfld->ipdtmpl[17]);
-            rec->LvlVal2 = rec->LvlVal1;
-
-         //  check for special fixed level types (1 through 10 or 101) and set the level values to 0
-         //  Reference: https://www.nco.ncep.noaa.gov/pmb/docs/grib2/grib2_doc/grib2_table4-5.shtml
+            rec->LvlVal2 = rec->LvlVal1;    
+           //  check for special fixed level types (1 through 10 or 101) and set the level values to 0
+           //  Reference: https://www.nco.ncep.noaa.gov/pmb/docs/grib2/grib2_doc/grib2_table4-5.shtml
          } else if( (rec->LvlTyp >= 1 && rec->LvlTyp <= 10) || rec->LvlTyp == 101 ) {
             rec->LvlVal1 = 0;
             rec->LvlVal2 = 0;
@@ -790,6 +802,11 @@ void MetGrib2DataFile::read_grib2_record_list() {
             rec->StatType = gfld->ipdtmpl[23];
          }
 
+         //  percentile value for templates 6 and 10
+         if( 6 == gfld->ipdtnum || 10 == gfld->ipdtnum ){
+            rec->PercVal = gfld->ipdtmpl[15];
+         }
+
          //  depending on the template number, determine the reference times
          if( 8 <= gfld->ipdtnum && 12 >= gfld->ipdtnum ){
 

src/libcode/vx_data2d_grib2/data2d_grib2.h

@@ -61,6 +61,7 @@ typedef struct {
    int EnsNumber;
    int DerType;
    int StatType;
+   int PercVal;
    IntArray IPDTmpl;
 } Grib2Record;
 

src/libcode/vx_data2d_grib2/var_info_grib2.cc

@@ -108,6 +108,7 @@ void VarInfoGrib2::assign(const VarInfoGrib2 &v) {
    EnsType    = v.EnsType;
    DerType    = v.DerType;
    StatType   = v.StatType;
+   PercVal    = v.PercVal;
 
    IPDTmplIndex = v.IPDTmplIndex;
    IPDTmplVal   = v.IPDTmplVal;
@@ -134,6 +135,7 @@ void VarInfoGrib2::clear() {
    EnsType    = bad_data_int;
    DerType    = bad_data_int;
    StatType   = bad_data_int;
+   PercVal    = bad_data_int;
 
    IPDTmplIndex.clear();
    IPDTmplVal.clear();
@@ -157,7 +159,8 @@ void VarInfoGrib2::dump(ostream &out) const {
        << "  Process    = " << Process    << "\n"
        << "  EnsType    = " << EnsType    << "\n"
        << "  DerType    = " << DerType    << "\n"
-       << "  StatType   = " << StatType   << "\n";
+       << "  StatType   = " << StatType   << "\n"
+       << "  PercVal    = " << PercVal    << "\n";
    out << "  IPDTmplIndex:\n";
    IPDTmplIndex.dump(out);
    out << "  IPDTmplVal:\n";
@@ -245,6 +248,13 @@ void VarInfoGrib2::set_stat_type(int v) {
 
 ///////////////////////////////////////////////////////////////////////////////
 
+void VarInfoGrib2::set_perc_val(int v) {
+   PercVal = v;
+   return;
+}
+
+///////////////////////////////////////////////////////////////////////////////
+
 void VarInfoGrib2::set_ipdtmpl_index(const IntArray &v) {
    IPDTmplIndex = v;
    return;
@@ -280,6 +290,7 @@ void VarInfoGrib2::set_dict(Dictionary & dict) {
    EnsType                 = dict.lookup_int   (conf_key_GRIB2_ens_type,  false);
    DerType                 = dict.lookup_int   (conf_key_GRIB2_der_type,  false);
    StatType                = dict.lookup_int   (conf_key_GRIB2_stat_type, false);
+   PercVal                 = dict.lookup_int   (conf_key_GRIB2_perc_val,  false);
 
    IPDTmplIndex = dict.lookup_int_array(conf_key_GRIB2_ipdtmpl_index, false);
    IPDTmplVal   = dict.lookup_int_array(conf_key_GRIB2_ipdtmpl_val,   false);

src/libcode/vx_data2d_grib2/var_info_grib2.h

@@ -56,6 +56,7 @@ class VarInfoGrib2 : public VarInfo
       int EnsType;    // Type of Ensemble Forecast (Table 4.6)
       int DerType;    // Derived Forecast (Table 4.7)
       int StatType;   // Statistical Processing Type (Table 4.10)
+      int PercVal;    // Percentile Value (Octet 35 for Templates 4.6 and 4.10)
 
       IntArray IPDTmplIndex; // Index into the GRIB2 ipdtmpl array
       IntArray IPDTmplVal;   // Corresponding GRIB2 ipdtmpl value
@@ -88,6 +89,7 @@ class VarInfoGrib2 : public VarInfo
       int         ens_type()    const;
       int         der_type()    const;
       int         stat_type()   const;
+      int         perc_val()    const;
 
       int         n_ipdtmpl()        const;
       int         ipdtmpl_index(int) const;
@@ -110,6 +112,7 @@ class VarInfoGrib2 : public VarInfo
       void set_ens_type(int);
       void set_der_type(int);
       void set_stat_type(int);
+      void set_perc_val(int);
       void set_ipdtmpl_index(const IntArray &);
       void set_ipdtmpl_val(const IntArray &);
 
@@ -143,6 +146,7 @@ inline int         VarInfoGrib2::process()    const { return(Process);      }
 inline int         VarInfoGrib2::ens_type()   const { return(EnsType);      }
 inline int         VarInfoGrib2::der_type()   const { return(DerType);      }
 inline int         VarInfoGrib2::stat_type()  const { return(StatType);     }
+inline int         VarInfoGrib2::perc_val()   const { return(PercVal);      }
 inline int         VarInfoGrib2::n_ipdtmpl()  const {
                                     return(IPDTmplIndex.n()); }
 inline int         VarInfoGrib2::ipdtmpl_index(int i) const {