Bugfix #2833 develop azimuth (#2840)

* Per #2833, fix n-1 bug when defining the azimuth delta for range/azimuth grids.

* Per #2833, when definng TcrmwData:range_max_km, divide by n_range - 1 since the range values start at 0.

* Per #2833, remove max_range_km from the TC-RMW config file. Set the default rmw_scale to NA so that its not used by default. And update the documentation. Still actually need to make the logic of the code work as it should.

* Per #2833, update tc_rmw to define the range as either a function of rmw or using explicit spacing in km.

* Per #2833, update the TCRMW Config files to remove the max_range_km entry, and update the unit test for one call to use RMW ranges and the other to use ranges defined in kilometers.

* Per #2833, just correct code comments.

* Per #2833, divide by n - 1 when computing the range delta, rather than n.

* Per #2833, correct the handling of the maximum range in the tc-rmw tool. For fixed delta km, need to define the max range when setting up the grid at the beginning.

---------

Co-authored-by: MET Tools Test Account <met_test@seneca.rap.ucar.edu>

data/config/TCRMWConfig_default

@@ -98,9 +98,8 @@ regrid = {
 //
 n_range        = 100;
 n_azimuth      = 180;
-max_range_km   = 1000.0;
 delta_range_km = 10.0;
-rmw_scale      = 0.2;
+rmw_scale      = NA;
 
 //
 // Optionally convert u/v winds to tangential/radial winds

docs/Users_Guide/tc-rmw.rst

@@ -101,27 +101,19 @@ The **n_azimuth** parameter is the number of equally spaced azimuth intervals in
 
 _______________________
 
-.. code-block:: none
-
-  max_range_km = 100.0;
-
-The **max_range_km** parameter specifies the maximum range of the range-azimuth grid, in kilometers. If this parameter is specified and not **rmw_scale**, the radial grid spacing will be **max_range_km / n_range**.
-
-_______________________
-
 .. code-block:: none
 
   delta_range_km = 10.0;
 
-The **delta_range_km** parameter specifies the spacing of the range rings, in kilometers.
+The **delta_range_km** parameter specifies the spacing of the range rings, in kilometers. The range values start with 0 km and extend out to **n_range - 1** times this delta spacing.
 
 _______________________
 
 .. code-block:: none
 
-  rmw_scale = 0.2;
+  rmw_scale = NA;
 
-The **rmw_scale** parameter overrides the **max_range_km** parameter. When this is set the radial grid spacing will be **rmw_scale** in units of the RMW, which varies along the storm track.
+If changed from its default value of **NA**, the **rmw_scale** parameter overrides the **delta_range_km** parameter. The radial grid spacing is defined using **rmw_scale** in units of the RMW, which varies along the storm track. For example, setting **rmw_scale** to 0.2 would define the delta range spacing as 20% of the radius of maximum winds around each point. Note that RMW is defined in nautical miles but is converted to kilometers for this computation. 
 
 _______________________
 

internal/test_unit/config/TCRMWConfig_gonzalo

@@ -99,9 +99,8 @@ regrid = {
 //
 n_range        = 50;
 n_azimuth      = 90;
-max_range_km   = 1000.0;
 delta_range_km = 10.0;
-rmw_scale      = 0.2;
+rmw_scale      = NA;
 
 //
 // Optionally convert u/v winds to tangential/radial winds

internal/test_unit/config/TCRMWConfig_pressure_lev_out

@@ -99,8 +99,7 @@ regrid = {
 //
 n_range        = 100;
 n_azimuth      = 180;
-max_range_km   = 1000.0;
-delta_range_km = 10.0;
+delta_range_km = NA;
 rmw_scale      = 0.2;
 
 //

src/basic/vx_config/config_constants.h

@@ -1213,7 +1213,6 @@ static const char conf_key_nc_pairs_grid[]                   = "nc_pairs_grid";
 
 static const char conf_key_n_range[]     = "n_range";
 static const char conf_key_n_azimuth[]   = "n_azimuth";
-static const char conf_key_max_range[]   = "max_range_km";
 static const char conf_key_delta_range[] = "delta_range_km";
 static const char conf_key_rmw_scale[]   = "rmw_scale";
 static const char conf_key_compute_tangential_and_radial_winds[] = "compute_tangential_and_radial_winds";

src/libcode/vx_grid/tcrmw_grid.cc

@@ -213,8 +213,9 @@ RLLD.rot_lat_ll = 90.0 - range_max_deg;
 RLLD.rot_lon_ll =  0.0;
 
 RLLD.delta_rot_lat = range_max_deg/(Range_n - 1);
-// RLLD.delta_rot_lon = 360.0/Azimuth_n;
-RLLD.delta_rot_lon = 360.0/(Azimuth_n - 1);
+
+// MET #2833 divide by n rather than n-1 for the azimuth increment
+RLLD.delta_rot_lon = 360.0/Azimuth_n;
 
 RLLD.Nlat = Range_n;
 RLLD.Nlon = Azimuth_n;

src/libcode/vx_grid/tcrmw_grid.h

@@ -69,7 +69,7 @@ class TcrmwGrid : public RotatedLatLonGrid {
       int azimuth_n () const;
 
       double range_max_km      () const;
-      double range_delta_km    () const;   //  Range_Max_km/Range_n
+      double range_delta_km    () const;   //  Range_max_km/(Range_n - 1)
 
       double azimuth_delta_deg () const;   //  360.0/Azimuth_n
 
@@ -117,7 +117,7 @@ inline int TcrmwGrid::azimuth_n () const { return ( Azimuth_n ); }
 
 inline double TcrmwGrid::range_max_km () const { return ( Range_max_km ); }
 
-inline double TcrmwGrid::range_delta_km () const { return ( Range_max_km/Range_n ); }
+inline double TcrmwGrid::range_delta_km () const { return ( Range_max_km/(Range_n - 1) ); }
 
 inline double TcrmwGrid::azimuth_delta_deg () const { return ( 360.0/Azimuth_n ); }
 

src/tools/tc_utils/tc_diag/tc_diag.cc

@@ -17,6 +17,7 @@
 //   000    09/27/22  Halley Gotway New
 //   001    08/17/23  Halley Gotway MET #2609 handle missing data
 //   002    10/24/23  Halley Gotway MET #2550 enhance diagnostics
+//   003    03/11/24  Halley Gotway MET #2833 range/azimuth grid
 //
 ////////////////////////////////////////////////////////////////////////
 
@@ -2203,7 +2204,9 @@ void TmpFileInfo::setup_nc_file(const DomainInfo &di,
    // Set grid center
    d.lat_center   =      pnt_ptr->lat();
    d.lon_center   = -1.0*pnt_ptr->lon(); // degrees east to west
-   d.range_max_km = di.delta_range_km * d.range_n;
+
+   // MET #2833 multiply by n-1 since the ranges begin at 0 km
+   d.range_max_km = di.delta_range_km * (d.range_n - 1);
 
    // Instantiate the grid
    grid_out.set(d);

src/tools/tc_utils/tc_rmw/tc_rmw.cc

@@ -19,6 +19,7 @@
 //   002   07/06/22  Howard Soh      METplus-Internal #19 Rename main to met_main
 //   003   09/28/22  Prestopnik      MET #2227 Remove namspace std and netCDF from header files
 //   004   04/26/23  Halley Gotway   MET #2523 Reorder NetCDF dimensions
+//   005   03/11/24  Halley Gotway   MET #2833 range/azimuth grid
 //
 ////////////////////////////////////////////////////////////////////////
 
@@ -554,7 +555,13 @@ void setup_grid() {
     grid_data.name = "TCRMW";
     grid_data.range_n = conf_info.n_range;
     grid_data.azimuth_n = conf_info.n_azimuth;
-    grid_data.range_max_km = conf_info.max_range_km;
+
+    // Define the maximum range in km based on the fixed increment 
+    if(is_bad_data(conf_info.rmw_scale)) {
+        grid_data.range_max_km =
+            conf_info.delta_range_km *
+            (conf_info.n_range - 1);
+    }
 
     tcrmw_grid.set_from_data(grid_data);
     grid.set(grid_data);
@@ -593,8 +600,9 @@ void setup_nc_file() {
         lead_time_str_var, lead_time_sec_var);
 
     // Define range and azimuth dimensions
-    def_tc_range_azimuth(nc_out, range_dim, azimuth_dim, tcrmw_grid,
-        conf_info.rmw_scale);
+    def_tc_range_azimuth(nc_out,
+        range_dim, azimuth_dim,
+        tcrmw_grid, conf_info.rmw_scale);
 
     // Define latitude and longitude arrays
     def_tc_lat_lon(nc_out,
@@ -606,11 +614,11 @@ void setup_nc_file() {
         // Get VarInfo
         data_info = conf_info.data_info[i_var];
         mlog << Debug(4) << "Processing field: " << data_info->magic_str() << "\n";
-	     string fname = data_info->name_attr();
+        string fname = data_info->name_attr();
         variable_levels[fname].push_back(data_info->level_attr());
         variable_long_names[fname] = data_info->long_name_attr();
         variable_units[fname] = data_info->units_attr();
-	     wind_converter.update_input(fname, data_info->units_attr());
+        wind_converter.update_input(fname, data_info->units_attr());
     }
 
     // Define pressure levels
@@ -644,7 +652,7 @@ void compute_lat_lon(TcrmwGrid& tcrmw_grid,
                 ia * tcrmw_grid.azimuth_delta_deg(),
                 lat, lon);
             lat_arr[i] = lat;
-            lon_arr[i] = - lon;
+            lon_arr[i] = -lon;
         }
     }
 }
@@ -688,9 +696,15 @@ void process_fields(const TrackInfoArray& tracks) {
         grid_data.lat_center = point.lat();
         grid_data.lon_center = -1.0*point.lon(); // internal sign change
 
-        // RMW is same as mrd()
-        grid_data.range_max_km = conf_info.rmw_scale *
-            point.mrd() * tc_km_per_nautical_miles * conf_info.n_range;
+        // Define the maximum range in km relative to the radius of maximum winds
+        if(!is_bad_data(conf_info.rmw_scale)) {
+            grid_data.range_max_km =
+                conf_info.rmw_scale *
+                point.mrd() * tc_km_per_nautical_miles *
+                (conf_info.n_range - 1);
+        }
+
+        // Re-define the range/azimuth grid
         tcrmw_grid.clear();
         tcrmw_grid.set_from_data(grid_data);
         grid.clear();
@@ -713,7 +727,7 @@ void process_fields(const TrackInfoArray& tracks) {
 
         for(int i_var = 0; i_var < conf_info.get_n_data(); i_var++) {
 
-            // Update the variable info with the valid time of the track point
+            // Update with the valid time of the track point
             data_info = conf_info.data_info[i_var];
 
             string sname = data_info->name_attr().string();
@@ -731,24 +745,27 @@ void process_fields(const TrackInfoArray& tracks) {
             mlog << Debug(4) << "data_max:" << data_max << "\n";
 
             // Regrid data
-            data_dp = met_regrid(data_dp, latlon_arr, grid, data_info->regrid());
+            data_dp = met_regrid(data_dp, latlon_arr, grid,
+                                 data_info->regrid());
             data_dp.data_range(data_min, data_max);
             mlog << Debug(4) << "data_min:" << data_min << "\n";
             mlog << Debug(4) << "data_max:" << data_max << "\n";
 
-	         // if this is "U", setup everything for matching "V" and compute the radial/tangential
-	         if(wind_converter.compute_winds_if_input_is_u(i_point, sname, slevel, valid_time, data_files, ftype,
-                   latlon_arr, lat_arr, lon_arr, grid, data_dp, tcrmw_grid)) {
-                write_tc_pressure_level_data(nc_out, tcrmw_grid,
-                    pressure_level_indices, data_info->level_attr(), i_point,
-					     data_3d_vars[conf_info.radial_velocity_field_name.string()],
-					     wind_converter.get_wind_r_arr());
-	             write_tc_pressure_level_data(nc_out, tcrmw_grid,
-					     pressure_level_indices, data_info->level_attr(), i_point,
-					     data_3d_vars[conf_info.tangential_velocity_field_name.string()],
-					     wind_converter.get_wind_t_arr());
+            // If this is "U", setup everything for matching "V"
+            // and compute the radial/tangential winds
+           if(wind_converter.compute_winds_if_input_is_u(
+                  i_point, sname, slevel, valid_time, data_files, ftype,
+                  latlon_arr, lat_arr, lon_arr, grid, data_dp, tcrmw_grid)) {
+               write_tc_pressure_level_data(nc_out, tcrmw_grid,
+                   pressure_level_indices, data_info->level_attr(), i_point,
+                   data_3d_vars[conf_info.radial_velocity_field_name.string()],
+                   wind_converter.get_wind_r_arr());
+               write_tc_pressure_level_data(nc_out, tcrmw_grid,
+                   pressure_level_indices, data_info->level_attr(), i_point,
+                   data_3d_vars[conf_info.tangential_velocity_field_name.string()],
+                   wind_converter.get_wind_t_arr());
             }
-	    
+
             // Write data
             if(variable_levels[data_info->name_attr()].size() > 1) {
                 write_tc_pressure_level_data(nc_out, tcrmw_grid,

src/tools/tc_utils/tc_rmw/tc_rmw.h

@@ -146,10 +146,6 @@ static Grid      grid;
 static double* lat_arr;
 static double* lon_arr;
 
-// Wind arrays
-/* static double* wind_r_arr; */
-/* static double* wind_t_arr; */
-
 ////////////////////////////////////////////////////////////////////////
 
 #endif  //  __TC_RMW_H__

src/tools/tc_utils/tc_rmw/tc_rmw_conf_info.cc

@@ -72,7 +72,6 @@ void TCRMWConfInfo::clear() {
 
     n_range        = bad_data_int;
     n_azimuth      = bad_data_int;
-    max_range_km   = bad_data_double;
     delta_range_km = bad_data_double;
     rmw_scale      = bad_data_double;
 
@@ -175,15 +174,21 @@ void TCRMWConfInfo::process_config(GrdFileType ftype) {
     // Conf: n_azimuth
     n_azimuth = Conf.lookup_int(conf_key_n_azimuth);
 
-    // Conf: max_range
-    max_range_km = Conf.lookup_double(conf_key_max_range);
-
     // Conf: delta_range
     delta_range_km = Conf.lookup_double(conf_key_delta_range);
 
     // Conf: rmw_scale
     rmw_scale = Conf.lookup_double(conf_key_rmw_scale);
 
+    // Error check
+    if(is_bad_data(delta_range_km) && is_bad_data(rmw_scale)) {
+        mlog << Error << "\nTCRMWConfInfo::process_config() -> "
+             << "the \"" << conf_key_delta_range << "\" and \""
+             << conf_key_rmw_scale << "\" configuration options "
+             << "cannot both be set to bad data.\n\n";
+        exit(1);
+    }
+
     compute_tangential_and_radial_winds = Conf.lookup_bool(conf_key_compute_tangential_and_radial_winds);
     u_wind_field_name = Conf.lookup_string(conf_key_u_wind_field_name);
     v_wind_field_name = Conf.lookup_string(conf_key_v_wind_field_name);
@@ -192,7 +197,6 @@ void TCRMWConfInfo::process_config(GrdFileType ftype) {
     tangential_velocity_long_field_name = Conf.lookup_string(conf_key_tangential_velocity_long_field_name);
     radial_velocity_long_field_name = Conf.lookup_string(conf_key_radial_velocity_long_field_name);
 
-
     // Conf: data.field
     fdict = Conf.lookup_array(conf_key_data_field);
 

src/tools/tc_utils/tc_rmw/tc_rmw_conf_info.h

@@ -50,18 +50,17 @@ class TCRMWConfInfo {
         // Range/Azimuth information
         int    n_range;
         int    n_azimuth;
-        double max_range_km;
         double delta_range_km;
         double rmw_scale;
 
-	// Wind conversion information
-	bool compute_tangential_and_radial_winds;
-	ConcatString u_wind_field_name;
-	ConcatString v_wind_field_name;
-	ConcatString tangential_velocity_field_name;
-	ConcatString radial_velocity_field_name;
-	ConcatString tangential_velocity_long_field_name;
-	ConcatString radial_velocity_long_field_name;
+        // Wind conversion information
+        bool compute_tangential_and_radial_winds;
+        ConcatString u_wind_field_name;
+        ConcatString v_wind_field_name;
+        ConcatString tangential_velocity_field_name;
+        ConcatString radial_velocity_field_name;
+        ConcatString tangential_velocity_long_field_name;
+        ConcatString radial_velocity_long_field_name;
 
         // Variable information
         VarInfo** data_info;