Merge pull request #1297 from NCAR/feature_1293_regrid

Pull request to merge feature_1293_regrid into develop

met/data/config/README

@@ -410,13 +410,19 @@ obtype = "ANALYS";
 //
 //   - The "gaussian_dx" entry specifies a delta distance for Gaussian
 //     smoothing. The default is 81.271. Ignored if not Gaussian method.
-// 
+//
 //   - The "gaussian_radius" entry defines the radius of influence for Gaussian
 //     smoothing. The default is 120. Ignored if not Gaussian method.
 //
-//   - The "gaussian_dx" and "gaussian_radius" settings must be in the same units,
-//     such as kilometers or degress.  Their ratio (sigma = gaussian_radius /
-//     gaussian_dx) determines the Guassian weighting function.
+//   - The "gaussian_dx" and "gaussian_radius" settings must be in the same
+//     units, such as kilometers or degress.  Their ratio
+//     (sigma = gaussian_radius / gaussian_dx) determines the Guassian weighting
+//     function.
+//
+//   - The "convert", "censor_thresh", and "censor_val" entries are described
+//     below.  When specified, these operations are applied to the output of the
+//     regridding step.  The conversion operation is applied first, followed by
+//     the censoring operation.
 //
 regrid = {
    to_grid         = NONE;
@@ -426,6 +432,9 @@ regrid = {
    shape           = SQUARE;
    gaussian_dx     = 81.271;
    gaussian_radius = 120;
+   convert(x)      = x;
+   censor_thresh   = [];
+   censor_val      = [];
 }
 
 //
@@ -921,7 +930,7 @@ climo_stdev = {
 }
 
 //
-// The "climo_cdf" dictionary specifies how the the climatological mean 
+// The "climo_cdf" dictionary specifies how the the climatological mean
 // ("climo_mean") and standard deviation ("climo_stdev") data are used to
 // evaluate model performance relative to where the observation value falls
 // within the climatological distribution. This dictionary consists of 3
@@ -941,7 +950,7 @@ climo_stdev = {
 // climatological bins. The array must begin with 0.0 and end with 1.0.
 // For example:
 //   cdf_bins = [ 0.0, 0.10, 0.25, 0.75, 0.90, 1.0 ];
-// 
+//
 // When "cdf_bins" is set to an integer, it defines the number of bins to be
 // used. The "center_bins" flag indicates whether or not the bins should be
 // centered on 0.5. An odd number of bins can be centered or uncentered while
@@ -980,7 +989,7 @@ climo_cdf = {
 //
 // When "climo_mean" and "climo_stdev" are both set to non-probability fields,
 // the MET tools use the mean, standard deviation, and observation event
-// threshold to derive a normal approximation of the climatological 
+// threshold to derive a normal approximation of the climatological
 // probabilities. Those derived probability values are used to compute BSS.
 //
 

met/src/basic/vx_config/config_constants.h

@@ -11,10 +11,12 @@
 #define __CONFIG_CONSTANTS_H__
 
 #include "vx_util.h"
+
 #include "GridTemplate.h"
 #include "int_array.h"
 #include "gsl_randist.h"
 #include "config_gaussian.h"
+#include "config_funcs.h"
 
 ////////////////////////////////////////////////////////////////////////
 
@@ -271,15 +273,20 @@ struct RegridInfo {
                             // or explicit grid definition.
    InterpMthd   method;     // Regridding method
    int          width;      // Regridding width
-   GaussianInfo gaussian;  // Gaussian smoothing
+   GaussianInfo gaussian;   // Gaussian smoothing
    GridTemplateFactory::GridTemplates shape; // Interpolation shape
    RegridInfo();
 
-   void *       hook;       // not allocated
+   // Process the regridded data
+   UserFunc_1Arg convert_fx;    // Conversion function
+   ThreshArray   censor_thresh; // Censoring thesholds
+   NumArray      censor_val;    // and replacement values
+
+   void *       hook;           // not allocated
 
    void         clear();
-   void         validate(); // ensure that width and method are accordant
-   void         validate_point();   // ensure that width and method are accordant
+   void         validate();        // ensure that width and method are accordant
+   void         validate_point();  // ensure that width and method are accordant
 };
 
 ////////////////////////////////////////////////////////////////////////

met/src/basic/vx_config/config_util.cc

@@ -138,6 +138,9 @@ void RegridInfo::clear() {
    width = bad_data_int;
    gaussian.clear();
    shape = GridTemplateFactory::GridTemplate_None;
+   convert_fx.clear(); 
+   censor_thresh.clear();
+   censor_val.clear();
 }
 
 ///////////////////////////////////////////////////////////////////////////////
@@ -209,13 +212,23 @@ void RegridInfo::validate() {
    // Check the Gaussian filter
    if(method == InterpMthd_MaxGauss) {
       if(gaussian.radius < gaussian.dx) {
-         mlog << Error << "\n"
+         mlog << Error << "\nRegridInfo::validate() -> "
               << "The radius of influence (" << gaussian.radius
               << ") is less than the delta distance (" << gaussian.dx
               << ") for regridding method \"" << interpmthd_to_string(method) << "\".\n\n";
          exit(1);
       }
    }
+
+   // Check for equal number of censor thresholds and values
+   if(censor_thresh.n() != censor_val.n()) {
+      mlog << Error << "\nRegridInfo::validate() -> "
+           << "The number of censor thresholds in \""
+           << conf_key_censor_thresh << "\" (" << censor_thresh.n()
+           << ") must match the number of replacement values in \""
+           << conf_key_censor_val << "\" (" << censor_val.n() << ").\n\n";
+      exit(1);
+   }
 
 }
 
@@ -1169,6 +1182,16 @@ RegridInfo parse_conf_regrid(Dictionary *dict, bool error_out) {
    info.gaussian.trunc_factor = (is_bad_data(conf_value) ? default_trunc_factor : conf_value);
    if (info.method == InterpMthd_Gaussian) info.gaussian.compute();
 
+   // Conf: convert
+   info.convert_fx.set(regrid_dict->lookup(conf_key_convert));
+
+   // Conf: censor_thresh
+   info.censor_thresh = regrid_dict->lookup_thresh_array(conf_key_censor_thresh, false);
+
+   // Conf: censor_val
+   info.censor_val = regrid_dict->lookup_num_array(conf_key_censor_val, false);
+
+   // Validate the settings
    info.validate();
 
    return(info);

met/src/basic/vx_util/data_plane.cc

@@ -282,6 +282,21 @@ void DataPlane::threshold(const SingleThresh &st) {
 
 ///////////////////////////////////////////////////////////////////////////////
 
+void DataPlane::convert(const UserFunc_1Arg &convert_fx) {
+
+   if(!convert_fx.is_set()) return;
+
+   mlog << Debug(3) << "Applying conversion function.\n";
+
+   for(int i=0; i<Nxy; i++) {
+      if(!is_bad_data(buf()[i])) buf()[i] = convert_fx(buf()[i]);
+   }
+
+   return;
+}
+
+///////////////////////////////////////////////////////////////////////////////
+
 void DataPlane::censor(const ThreshArray &censor_thresh,
                        const NumArray &censor_val) {
    int i, j, count;

met/src/basic/vx_util/data_plane.h

@@ -20,6 +20,7 @@
 #include "thresh_array.h"
 #include "num_array.h"
 #include "two_d_array.h"
+#include "config_funcs.h"
 
 #include "vx_cal.h"
 
@@ -100,6 +101,7 @@ class DataPlane {
          //
 
       void threshold(const SingleThresh &);
+      void convert  (const UserFunc_1Arg &);
       void censor   (const ThreshArray &, const NumArray &);
 
       void replace_bad_data(const double value);

met/src/libcode/vx_data2d/data_class.cc

@@ -285,18 +285,7 @@ if ( ShiftRight != 0 )  dp.shift_right(ShiftRight);
    // Apply conversion logic.
    //
 
-if ( vinfo->ConvertFx.is_set() )  {
-
-   mlog << Debug(3) << "Applying conversion function.\n";
-
-   int Nxy = dp.nx()*dp.ny();
-
-   for (int j=0; j<Nxy; ++j)  {
-      if ( ! is_bad_data(dp.buf()[j]) )  {
-         dp.buf()[j] = vinfo->ConvertFx(dp.buf()[j]);
-      }
-   }
-}
+dp.convert(vinfo->ConvertFx);
 
    //
    // Apply censor logic.

met/src/libcode/vx_regrid/vx_regrid.cc

@@ -68,6 +68,18 @@ switch ( info.method )  {
 
 }   //  switch info.method
 
+   //
+   //  apply convert logic
+   //
+
+out.convert(info.convert_fx);
+
+   //
+   //  apply censor logic
+   //
+
+out.censor(info.censor_thresh, info.censor_val);
+
    //
    //  done
    //

met/src/libcode/vx_statistics/pair_base.cc

@@ -801,7 +801,7 @@ double PairBase::process_obs(VarInfo *vinfo, double v) {
 
 ////////////////////////////////////////////////////////////////////////
 //
-// Begin miscellanous utility functions
+// Begin miscellaneous utility functions
 //
 ////////////////////////////////////////////////////////////////////////
 

met/src/tools/other/mode_graphics/Makefile.am

@@ -35,6 +35,7 @@ plot_mode_field_LDADD = -lvx_config \
 	-lvx_render \
 	-lvx_color \
 	-lvx_util \
+	-lvx_config \
 	-lvx_cal \
 	-lvx_log \
 	-lvx_math \

met/src/tools/other/regrid_data_plane/regrid_data_plane.cc

@@ -22,6 +22,7 @@
 //   002    06-25-17  Howard Soh     Support GOES-16
 //   003    09-24-17  Howard Soh     Support Gaussian filtering
 //   004    01-28-20  Howard Soh     Moved GOES-16/17 to point2grib
+//   005    04-09-20  Halley Gotway  Add convert and censor options.
 //
 ////////////////////////////////////////////////////////////////////////
 
@@ -76,7 +77,6 @@ static NcFile *nc_out  = (NcFile *) 0;
 static NcDim  lat_dim ;
 static NcDim  lon_dim ;
 
-
 ////////////////////////////////////////////////////////////////////////
 
 static void process_command_line(int, char **);
@@ -93,6 +93,8 @@ static void set_gaussian_dx(const StringArray &);
 static void set_gaussian_radius(const StringArray &);
 static void set_width(const StringArray &);
 static void set_vld_thresh(const StringArray &);
+static void set_convert(const StringArray &);
+static void set_censor(const StringArray &);
 static void set_name(const StringArray &);
 static void set_logfile(const StringArray &);
 static void set_verbosity(const StringArray &);
@@ -150,6 +152,8 @@ void process_command_line(int argc, char **argv) {
    cline.add(set_gaussian_radius, "-gaussian_radius", 1);
    cline.add(set_gaussian_dx,     "-gaussian_dx",      1);
    cline.add(set_vld_thresh, "-vld_thresh", 1);
+   cline.add(set_convert,    "-convert",    1);
+   cline.add(set_censor,     "-censor",     2);
    cline.add(set_name,       "-name",       1);
    cline.add(set_logfile,    "-log",        1);
    cline.add(set_verbosity,  "-v",          1);
@@ -464,6 +468,8 @@ void usage() {
         << "\t[-gaussian_radius n]\n"
         << "\t[-shape type]\n"
         << "\t[-vld_thresh n]\n"
+        << "\t[-convert string]\n"
+        << "\t[-censor thresh value]\n"
         << "\t[-name list]\n"
         << "\t[-log file]\n"
         << "\t[-v level]\n"
@@ -472,8 +478,8 @@ void usage() {
         << "\twhere\t\"input_filename\" is the gridded data file to be "
         << "read (required).\n"
 
-        << "\t\t\"to_grid\" defines the output grid as a named grid, the "
-        << "path to a gridded data file, or an explicit grid "
+        << "\t\t\"to_grid\" defines the output grid as a named grid, "
+        << "the path to a gridded data file, or an explicit grid "
         << "specification string (required).\n"
 
         << "\t\t\"output_filename\" is the output NetCDF file to be "
@@ -488,24 +494,35 @@ void usage() {
 
         << "\t\t\"-width n\" overrides the default regridding "
         << "width (" << RGInfo.width << ") (optional).\n"
-        << "\t\t\tThe width should be the ratio of dx between from_grid and to_grid for MAXGAUSS.\n"
-        << "\t\t\tFor example, width=" << nint(RGInfo.gaussian.dx / 3) << " if the from_grid is 3 km and to_grid is "
+        << "\t\t\tThe width should be the ratio of dx between "
+        << "from_grid and to_grid for MAXGAUSS.\n"
+        << "\t\t\tFor example, width=" << nint(RGInfo.gaussian.dx / 3)
+        << " if the from_grid is 3 km and to_grid is "
         << RGInfo.gaussian.dx << "km.\n"
 
-        << "\t\t\"-gaussian_dx n\" specifies a delta distance for Gaussian smoothing."
-        << " The default is " << RGInfo.gaussian.dx << ". Ignored if not Gaussian method (optional).\n"
+        << "\t\t\"-gaussian_dx n\" overrides the default a delta "
+        << "distance for Gaussian smoothing (" << RGInfo.gaussian.dx
+        << ") (optional).\n"
 
-        << "\t\t\"-gaussian_radius n\" specifies the radius of influence for Gaussian smoothing."
-        << " The default is " << RGInfo.gaussian.radius << ". Ignored if not Gaussian method (optional).\n"
+        << "\t\t\"-gaussian_radius n\" overrides the default radius of "
+        << "influence for Gaussian smoothing ("
+        << RGInfo.gaussian.radius << ") (optional).\n"
 
-        << "\t\t\"-shape type\" overrides the default interpolation shape ("
-        << gtf.enum2String(RGInfo.shape) << ") "
-        <<  "(optional).\n"
+        << "\t\t\"-shape type\" overrides the default interpolation "
+        << "shape (" << gtf.enum2String(RGInfo.shape) << ") (optional).\n"
 
         << "\t\t\"-vld_thresh n\" overrides the default required "
         << "ratio of valid data for regridding (" << RGInfo.vld_thresh
         << ") (optional).\n"
 
+        << "\t\t\"-convert string\" specifies a conversion for the "
+        << "regridded output (optional).\n"
+        << "\t\t\tFor example, -convert 'convert(x) = x - 273.15;'\n"
+
+        << "\t\t\"-censor thresh value\" specifies censoring logic for "
+        << "the regridded output as a threshold string and replacement "
+        << "value (optional).\n"
+
         << "\t\t\"-name list\" specifies a comma-separated list of "
         << "output variable names for each field specified (optional).\n"
 
@@ -515,7 +532,8 @@ void usage() {
         << "\t\t\"-v level\" overrides the default level of logging ("
         << mlog.verbosity_level() << ") (optional).\n"
 
-        << "\t\t\"-compress level\" overrides the compression level of NetCDF variable (optional).\n\n" << flush;
+        << "\t\t\"-compress level\" overrides the compression level of "
+        << "NetCDF variable (optional).\n\n" << flush;
 
    exit(1);
 }
@@ -571,6 +589,30 @@ void set_vld_thresh(const StringArray &a) {
 
 ////////////////////////////////////////////////////////////////////////
 
+void set_convert(const StringArray &a) {
+
+   // Can only be used once
+   if(RGInfo.convert_fx.is_set()) {
+      mlog << Error << "\nset_convert() -> "
+           << "-convert_x may only be used once!\n\n";
+      exit(1);
+   }
+
+   MetConfig config;
+   config.read_string(a[0].c_str());
+
+   RGInfo.convert_fx.set(config.lookup(conf_key_convert));
+}
+
+////////////////////////////////////////////////////////////////////////
+
+void set_censor(const StringArray &a) {
+   RGInfo.censor_thresh.add(a[0].c_str());
+   RGInfo.censor_val.add(atof(a[1].c_str()));
+}
+
+////////////////////////////////////////////////////////////////////////
+
 void set_name(const StringArray & a) {
    VarNameSA.add_css(a[0]);
 }