diff --git a/docs/Users_Guide/reformat_point.rst b/docs/Users_Guide/reformat_point.rst index fb0a7b2766..ae9691745f 100644 --- a/docs/Users_Guide/reformat_point.rst +++ b/docs/Users_Guide/reformat_point.rst @@ -1092,9 +1092,22 @@ Optional Arguments for point2grid 16. The **-compress level** option indicates the desired level of compression (deflate level) for NetCDF variables. The valid level is between 0 and 9. The value of "level" will override the default setting of 0 from the configuration file or the environment variable MET_NC_COMPRESS. Setting the compression level to 0 will make no compression for the NetCDF output. Lower number is for fast compression and higher number is for better compression. Only 4 interpolation methods are applied to the field variables; MIN/MAX/MEDIAN/UW_MEAN. The GAUSSIAN method is applied to the probability variable only. Unlike regrad_data_plane, MAX method is applied to the file variable and Gaussian method to the probability variable with the MAXGAUSS method. If the probability variable is not requested, MAXGAUSS method is the same as MAX method. - + For the GOES-16 and GOES-17 data, the computing lat/long is time consuming. The computed coordinate (lat/long) is saved to a temporary NetCDF file, as described in :numref:`Contributor's Guide Section %s `. The computing lat/long step can be skipped if the coordinate file is given through the environment variable MET_GEOSTATIONARY_DATA. The grid mapping to the target grid is saved to MET_TMP_DIR to save the execution time. Once this file is created, the MET_GEOSTATIONARY_DATA is ignored. The grid mapping file should be deleted manually in order to apply a new MET_GEOSTATIONARY_DATA environment variable or to re-generate the grid mapping file. An example of call point2grid to process GOES-16 AOD data is shown below: + +The grid name or the grid definition can be given with the -field option when the grid information is missing from the input NetCDF file for the latitude_longitude projection. The latitude and longitude variable names should be defined by the user, and the grid information from the set_attr_grid is ignored in this case + +.. code-block:: none + + point2grid \ + iceh.2018-01-03.c00.tlat_tlon.nc \ + G231 \ + point2grid_cice_to_G231.nc \ + -config Point2GridConfig_tlat_tlon \ + -field 'name="hi_d"; level="(0,*,*)"; set_attr_grid="latlon 1440 1080 -79.80672 60.28144 0.04 0.04";' \ + -v 1 + .. code-block:: none point2grid \ diff --git a/internal/test_unit/config/Point2GridConfig_tlat_tlon b/internal/test_unit/config/Point2GridConfig_tlat_tlon new file mode 100644 index 0000000000..66f2f854b4 --- /dev/null +++ b/internal/test_unit/config/Point2GridConfig_tlat_tlon @@ -0,0 +1,6 @@ +file_type = NETCDF_NCCF; + +var_name_map = [ + { key = "lat_vname"; val = "TLAT"; }, + { key = "lon_vname"; val = "TLON"; } +]; diff --git a/internal/test_unit/xml/unit_plot_data_plane.xml b/internal/test_unit/xml/unit_plot_data_plane.xml index 5c56a389d8..599bfd0742 100644 --- a/internal/test_unit/xml/unit_plot_data_plane.xml +++ b/internal/test_unit/xml/unit_plot_data_plane.xml @@ -680,4 +680,18 @@ + + &MET_BIN;/plot_data_plane + \ + &DATA_DIR_MODEL;/nccf/MITLL.ProxyEchoTopsCalibratedMosaic.20200831_235328_v_20200831_235328.nc \ + &OUTPUT_DIR;/plot_data_plane/EchoTops_set_attr_grid.ps \ + 'name="ProxyEchoTopsCalibratedMosaic"; level="(0,*,*)"; set_attr_grid="latlon 8008 4004 -90 -180 0.04 0.04";' \ + -title "Global Synthetic Weather Radar EchoTops" \ + -v 1 + + + &OUTPUT_DIR;/plot_data_plane/EchoTops_set_attr_grid.ps + + + diff --git a/internal/test_unit/xml/unit_point2grid.xml b/internal/test_unit/xml/unit_point2grid.xml index aa9a1d0410..840fbdc7ea 100644 --- a/internal/test_unit/xml/unit_point2grid.xml +++ b/internal/test_unit/xml/unit_point2grid.xml @@ -356,4 +356,23 @@ + + &MET_BIN;/point2grid + + MET_TMP_DIR &OUTPUT_DIR;/point2grid + + \ + &DATA_DIR_MODEL;/cice/iceh.2018-01-03.c00.tlat_tlon.nc \ + G231 \ + &OUTPUT_DIR;/point2grid/point2grid_cice_to_G231.nc \ + -config &CONFIG_DIR;/Point2GridConfig_tlat_tlon \ + -field 'name="hi_d"; level="(0,*,*)"; set_attr_grid="latlon 180 360 -80 -180 1 1";' \ + -v 1 + + + &OUTPUT_DIR;/point2grid/point2grid_cice_to_G231.nc + + + + diff --git a/src/basic/vx_log/vx_log.h b/src/basic/vx_log/vx_log.h index ce8a8132cb..0ae10c3811 100644 --- a/src/basic/vx_log/vx_log.h +++ b/src/basic/vx_log/vx_log.h @@ -26,6 +26,15 @@ //////////////////////////////////////////////////////////////////////// +inline double get_exe_duration(clock_t start_clock, clock_t end_clock) { + return ((double)(end_clock - start_clock)) / CLOCKS_PER_SEC; +} + +inline double get_exe_duration(clock_t start_clock) { + return get_exe_duration(start_clock, clock()); +} + +//////////////////////////////////////////////////////////////////////// #endif // __VX_LOG_H__ diff --git a/src/basic/vx_math/is_bad_data.h b/src/basic/vx_math/is_bad_data.h index d8d7b334c3..80af0916ef 100644 --- a/src/basic/vx_math/is_bad_data.h +++ b/src/basic/vx_math/is_bad_data.h @@ -48,7 +48,7 @@ inline int is_bad_data(float a) { } inline int is_bad_data(char a) { - return(a == bad_data_char); + return (a == bad_data_char); } inline int is_eq(double a, double b, double tol) { @@ -60,6 +60,30 @@ inline int is_eq(double a, double b) { return is_eq(a, b, default_tol); } +inline int is_eq(double a, int b) { + return is_eq(a, (double)b); +} + +inline int is_eq(int a, double b) { + return is_eq((double)a, b); +} + +inline int is_eq(double a, unixtime b) { + return is_eq(a, (double)b); +} + +inline int is_eq(unixtime a, double b) { + return is_eq((double)a, b); +} + +inline int is_eq(float a, float b) { + return is_eq((double)a, (double)b); +} + +template +inline int is_eq(T a, T b) { + return (a == b); +} //////////////////////////////////////////////////////////////////////// diff --git a/src/libcode/vx_data2d/data2d_utils.cc b/src/libcode/vx_data2d/data2d_utils.cc index 6284ca434b..6fc68a93cc 100644 --- a/src/libcode/vx_data2d/data2d_utils.cc +++ b/src/libcode/vx_data2d/data2d_utils.cc @@ -24,6 +24,54 @@ using namespace std; //////////////////////////////////////////////////////////////////////// +bool build_grid_by_grid_string(const char *grid_str, Grid &grid, + const char *caller_name, bool do_warning) { + bool status = false; + + if (nullptr != grid_str && m_strlen(grid_str) > 0) { + // Parse as a white-space separated string + + StringArray sa; + sa.parse_wsss(grid_str); + + // Search for a named grid + if (sa.n() == 1 && find_grid_by_name(sa[0].c_str(), grid)) { + status = true; + mlog << Debug(3) << "Use the grid named \"" << grid_str << "\".\n"; + } + // Parse grid definition + else if (sa.n() > 1 && parse_grid_def(sa, grid)) { + status = true; + mlog << Debug(3) << "Use the grid defined by string \"" + << grid_str << "\".\n"; + } + else if (do_warning) { + mlog << Warning << "\nbuild_grid_by_grid_string() by " << caller_name + << " unsupported " << conf_key_set_attr_grid + << " definition string (" << grid_str + << ")!\n\n"; + } + } + + return status; +} + +//////////////////////////////////////////////////////////////////////// + +bool build_grid_by_grid_string(const ConcatString &grid_str, Grid &grid, + const char *caller_name, bool do_warning) { + bool status = false; + + if(grid_str.nonempty()) { + status = build_grid_by_grid_string(grid_str.c_str(), grid, + caller_name, do_warning); + } + + return status; +} + +//////////////////////////////////////////////////////////////////////// + bool derive_wdir(const DataPlane &u2d, const DataPlane &v2d, DataPlane &wdir2d) { int x, y; diff --git a/src/libcode/vx_data2d/data2d_utils.h b/src/libcode/vx_data2d/data2d_utils.h index 317eb92f72..24299475d3 100644 --- a/src/libcode/vx_data2d/data2d_utils.h +++ b/src/libcode/vx_data2d/data2d_utils.h @@ -20,6 +20,14 @@ //////////////////////////////////////////////////////////////////////// +extern bool build_grid_by_grid_string(const char *attr_grid, Grid &grid, + const char *caller_name=nullptr, + bool do_warning=true); + +extern bool build_grid_by_grid_string(const ConcatString &attr_grid, Grid &grid, + const char *caller_name=nullptr, + bool do_warning=true); + extern bool derive_wdir(const DataPlane &u2d, const DataPlane &v2d, DataPlane &wdir); diff --git a/src/libcode/vx_data2d/data_class.cc b/src/libcode/vx_data2d/data_class.cc index 23ff40f93f..f364e34e42 100644 --- a/src/libcode/vx_data2d/data_class.cc +++ b/src/libcode/vx_data2d/data_class.cc @@ -257,13 +257,21 @@ mlog << Debug(3) << "Resetting grid definition from \"" // Make sure the grid dimensions do not change // - if ( raw_nx() != grid.nx() || raw_ny() != grid.ny() ) { + if ( raw_nx() <= 0 && raw_ny() <= 0 ) { - mlog << Error << "\nMet2dDataFile::set_grid() -> " + mlog << Warning << "\nMet2dDataFile::set_grid() -> " << "When resetting the grid definition to \"" << grid.serialize() << "\", the grid dimensions " - << "cannot change (" << grid.nx() << ", " << grid.ny() + << "are changed (" << grid.nx() << ", " << grid.ny() << ") != (" << raw_nx() << ", " << raw_ny() << ").\n\n"; + } + else if ( raw_nx() != grid.nx() || raw_ny() != grid.ny() ) { + + mlog << Error << "\nMet2dDataFile::set_grid() -> " + << "When resetting the grid definition to \"" + << grid.serialize() << "\", the grid dimensions " + << "cannot change to (" << grid.nx() << ", " << grid.ny() + << ") from (" << raw_nx() << ", " << raw_ny() << ").\n\n"; exit ( 1 ); diff --git a/src/libcode/vx_data2d/var_info.cc b/src/libcode/vx_data2d/var_info.cc index a1775a7275..2c92c6bf69 100644 --- a/src/libcode/vx_data2d/var_info.cc +++ b/src/libcode/vx_data2d/var_info.cc @@ -26,6 +26,7 @@ #include "vx_cal.h" #include "vx_math.h" #include "vx_log.h" +#include "data2d_utils.h" using namespace std; @@ -541,25 +542,7 @@ void VarInfo::set_dict(Dictionary &dict) { // Parse set_attr grid s = parse_set_attr_string(dict, conf_key_set_attr_grid); - if(s.nonempty()) { - - // Parse as a white-space separated string - StringArray sa; - sa.parse_wsss(s); - - // Search for a named grid - if(sa.n() == 1 && find_grid_by_name(sa[0].c_str(), SetAttrGrid)) { - } - // Parse grid definition - else if(sa.n() > 1 && parse_grid_def(sa, SetAttrGrid)) { - } - else { - mlog << Warning << "\nVarInfo::set_dict() -> " - << "unsupported " << conf_key_set_attr_grid - << " definition string (" << s - << ")!\n\n"; - } - } + build_grid_by_grid_string(s, SetAttrGrid, "VarInfo::set_dict(Dictionary &dict) ->"); // Parse set_attr times s = parse_set_attr_string(dict, conf_key_set_attr_init); diff --git a/src/libcode/vx_data2d_nc_cf/data2d_nc_cf.cc b/src/libcode/vx_data2d_nc_cf/data2d_nc_cf.cc index 8c4373a62c..cc77b7f17a 100644 --- a/src/libcode/vx_data2d_nc_cf/data2d_nc_cf.cc +++ b/src/libcode/vx_data2d_nc_cf/data2d_nc_cf.cc @@ -159,11 +159,14 @@ bool MetNcCFDataFile::data_plane(VarInfo &vinfo, DataPlane &plane) { // Not sure why we do this - NcVarInfo *data_var = (NcVarInfo *)nullptr; - VarInfoNcCF *vinfo_nc = (VarInfoNcCF *)&vinfo; + auto data_var = (NcVarInfo *)nullptr; + auto vinfo_nc = (VarInfoNcCF *)&vinfo; static const string method_name = "MetNcCFDataFile::data_plane(VarInfo &, DataPlane &) -> "; + Grid grid_attr = vinfo.grid_attr(); + _file->update_grid(grid_attr); + // Initialize the data plane plane.clear(); @@ -339,6 +342,9 @@ int MetNcCFDataFile::data_plane_array(VarInfo &vinfo, static const string method_name = "MetNcCFDataFile::data_plane_array(VarInfo &, DataPlaneArray &) -> "; + Grid grid_attr = vinfo.grid_attr(); + _file->update_grid(grid_attr); + // Initialize plane_array.clear(); diff --git a/src/libcode/vx_data2d_nc_cf/nc_cf_file.cc b/src/libcode/vx_data2d_nc_cf/nc_cf_file.cc index 3038965ab6..a739c5ae0d 100644 --- a/src/libcode/vx_data2d_nc_cf/nc_cf_file.cc +++ b/src/libcode/vx_data2d_nc_cf/nc_cf_file.cc @@ -57,8 +57,6 @@ static ConcatString y_dim_var_name; static double get_nc_var_att_double(const NcVar *nc_var, const char *att_name, bool is_required=true); -#define USE_BUFFER 1 - //////////////////////////////////////////////////////////////////////// @@ -136,6 +134,9 @@ void NcCfFile::close() _dims = (NcDim **)nullptr; } + grid_ready = false; + has_attr_grid = false; + _numDims = 0; _dimNames.clear(); @@ -213,8 +214,8 @@ bool NcCfFile::open(const char * filepath) // Pull out the variables int max_dim_count = 0; - NcVar *z_var = (NcVar *)nullptr; - NcVar *valid_time_var = (NcVar *)nullptr; + auto z_var = (NcVar *)nullptr; + auto valid_time_var = (NcVar *)nullptr; ConcatString att_value; StringArray varNames; @@ -311,7 +312,7 @@ bool NcCfFile::open(const char * filepath) // Parse the units for the time variable. ut = sec_per_unit = 0; if (get_var_units(valid_time_var, units)) { - if (units.length() == 0) { + if (units.empty()) { mlog << Warning << "\n" << method_name << "the \"time\" variable must contain a \"units\" attribute. " << "Using valid time of 0\n\n"; @@ -324,7 +325,7 @@ bool NcCfFile::open(const char * filepath) } NcVar bounds_time_var; - NcVar *nc_time_var = (NcVar *)nullptr; + auto nc_time_var = (NcVar *)nullptr; bool use_bounds_var = false; ConcatString bounds_var_name; nc_time_var = valid_time_var; @@ -341,10 +342,10 @@ bool NcCfFile::open(const char * filepath) if (bounds_att) delete bounds_att; // Determine the number of times present. - int n_times = (int) get_data_size(valid_time_var); + int n_times = get_data_size(valid_time_var); int tim_buf_size = n_times; if (use_bounds_var) tim_buf_size *= 2; - double *time_values = new double[tim_buf_size]; + auto time_values = new double[tim_buf_size]; if( get_nc_data(nc_time_var, time_values) ) { bool no_leap_year = get_att_no_leap_year(valid_time_var); @@ -403,7 +404,7 @@ bool NcCfFile::open(const char * filepath) // Parse the units for the time variable. if (get_var_units(&init_time_var, units)) { - if (units.length() == 0) { + if (units.empty()) { mlog << Warning << "\n" << method_name << "the \"forecast_reference_time\" variable must contain a \"units\" attribute.\n\n"; ut = sec_per_unit = 0; @@ -418,8 +419,8 @@ bool NcCfFile::open(const char * filepath) ut = sec_per_unit = 0; } - double time_value = get_nc_time(&init_time_var,(int)0); - InitTime = (unixtime)ut + sec_per_unit * time_value; + double time_value = get_nc_time(&init_time_var,0); + InitTime = ut + (unixtime)(sec_per_unit * time_value); } // Pull out the grid. This must be done after pulling out the dimension @@ -436,7 +437,8 @@ bool NcCfFile::open(const char * filepath) StringArray z_dims; StringArray t_dims; StringArray dimNames; - string var_x_dim_name, var_y_dim_name; + string var_x_dim_name; + string var_y_dim_name; if (IS_VALID_NC_P(_xDim)) var_x_dim_name = GET_NC_NAME_P(_xDim); if (IS_VALID_NC_P(_yDim)) var_y_dim_name = GET_NC_NAME_P(_yDim); for (int j=0; jx_slot) && (j != var->y_slot)) ) { - ++count; - if ( var == nullptr || ((j != var->x_slot) && (j != var->y_slot)) ) - { + if (has_attr_grid) { + mlog << Debug(3) << "\n" << method_name + << "star found in unknown slot (" << j << ") for " << GET_NC_NAME_P(v) << "\n\n"; + } + else { mlog << Error << "\n" << method_name << "star found in bad slot (" << j << ") for " << GET_NC_NAME_P(v) << "\n\n"; + exit(1); } } @@ -1011,13 +1041,22 @@ bool NcCfFile::getData(NcVar * v, const LongArray & a, DataPlane & plane) const } // check slots - additional logic to satisfy Fortify Null Dereference - int x_slot_tmp = 0; - int y_slot_tmp = 0; + int x_slot_tmp = 0; + int y_slot_tmp = 0; if (var == nullptr || var->x_slot < 0 || var->y_slot < 0) { - mlog << Error << "\n" << method_name - << "bad x|y|z slot\n\n"; - exit(1); + if (has_attr_grid) { + mlog << Warning << "\n" << method_name + << "bad x|y|z slot (" << var->x_slot << "|" << var->y_slot + << "|" << var->z_slot << "|" << var->t_slot <<")\n\n"; + x_slot_tmp = dim_count - 1; + y_slot_tmp = dim_count - 2; + } + else { + mlog << Error << "\n" << method_name + << "bad x|y|z slot\n\n"; + exit(1); + } } else { x_slot_tmp = var->x_slot; @@ -1116,7 +1155,7 @@ bool NcCfFile::getData(NcVar * v, const LongArray & a, DataPlane & plane) const // done mlog << Debug(6) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; return true; } @@ -1462,9 +1501,9 @@ void NcCfFile::read_netcdf_grid() !((_xDim && _yDim) || (x_dim_var_name.nonempty() && y_dim_var_name.nonempty()))) { - mlog << Error << "\nNcCfFile::read_netcdf_grid() -> " + mlog << Warning << "\nNcCfFile::read_netcdf_grid() -> " << "Couldn't figure out projection from information in netCDF file.\n\n"; - exit(1); + return; } return; @@ -1478,7 +1517,7 @@ void NcCfFile::read_netcdf_grid() Grid NcCfFile::build_grid_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, const long lat_counts, const long lon_counts) { Grid grid_ll; - bool swap_to_north; + bool swap_to_north = false; LatLonData data = get_data_from_lat_lon_vars(lat_var, lon_var, lat_counts, lon_counts, swap_to_north); @@ -1506,9 +1545,9 @@ void NcCfFile::get_grid_from_grid_mapping(const NcVarAtt *grid_mapping_att) bool status = get_att_value_chars(grid_mapping_att, mapping_name); if (!status) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Cannot extract grid mapping name from netCDF file.\n\n"; - exit(1); + return; } NcVar *grid_mapping_var = nullptr; @@ -1522,12 +1561,12 @@ void NcCfFile::get_grid_from_grid_mapping(const NcVarAtt *grid_mapping_att) } } /* endfor - i */ - if ((grid_mapping_var == 0) || (IS_INVALID_NC_P(grid_mapping_var))) + if ((nullptr == grid_mapping_var) || (IS_INVALID_NC_P(grid_mapping_var))) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Cannot extract grid mapping variable (" << mapping_name << ") from netCDF file.\n\n"; - exit(1); + return; } // Get the name of the grid mapping @@ -1536,9 +1575,10 @@ void NcCfFile::get_grid_from_grid_mapping(const NcVarAtt *grid_mapping_att) if (IS_INVALID_NC_P(grid_mapping_name_att)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Cannot get coordinate system name from netCDF file.\n\n"; - exit(1); + if (grid_mapping_name_att) delete grid_mapping_name_att; + return; } //string grid_mapping_name = grid_mapping_name_att->getValues(att->as_string(0); @@ -1606,10 +1646,10 @@ void NcCfFile::get_grid_from_grid_mapping(const NcVarAtt *grid_mapping_att) } else { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Unknown grid mapping name (" << grid_mapping_name << ") found in netCDF file.\n\n"; - exit(1); + return; } } @@ -1673,9 +1713,9 @@ void NcCfFile::get_grid_mapping_lambert_azimuthal_equal_area(const NcVar *grid_m x_coord_units_name == "meters") x_coord_to_m_cf = 1.0; else if (x_coord_units_name == "km") x_coord_to_m_cf = 1000.0; else { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "The X coordinates must be in meters or kilometers for MET.\n\n"; - exit(1); + return; } } } @@ -1698,9 +1738,9 @@ void NcCfFile::get_grid_mapping_lambert_azimuthal_equal_area(const NcVar *grid_m y_coord_units_name == "meters" ) y_coord_to_m_cf = 1.0; else if (y_coord_units_name == "km") y_coord_to_m_cf = 1000.0; else { - mlog << Error << "\n" << method_name << " -> " - << "The X coordinates must be in meters or kilometers for MET.\n\n"; - exit(1); + mlog << Warning << "\n" << method_name << " -> " + << "The Y coordinates must be in meters or kilometers for MET.\n\n"; + return; } } } @@ -1737,11 +1777,11 @@ void NcCfFile::get_grid_mapping_lambert_azimuthal_equal_area(const NcVar *grid_m double curr_delta = fabs(x_values[i] - x_values[i-1]); if (fabs(curr_delta - dx_m_a) > DELTA_TOLERANCE) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "MET can only process Lambert Azimuthal Equal Area files " << "where the delta along the x-axis is constant (" << curr_delta << " != " << dx_m_a << ")\n\n"; - exit(1); + return; } } @@ -1750,11 +1790,11 @@ void NcCfFile::get_grid_mapping_lambert_azimuthal_equal_area(const NcVar *grid_m double curr_delta = fabs(y_values[i] - y_values[i-1]); if (fabs(curr_delta - dy_m_a) > DELTA_TOLERANCE) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "MET can only process Lambert Azimuthal Equal Area files " << "where the delta along the y-axis is constant (" << curr_delta << " != " << dy_m_a << ")\n\n"; - exit(1); + return; } } @@ -1806,11 +1846,11 @@ void NcCfFile::get_grid_mapping_lambert_azimuthal_equal_area(const NcVar *grid_m if(!is_bad_data(false_easting) && !is_eq(false_easting, 0.0)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "MET cannot process Lambert Azimuthal Equal Area files " << "with non-zero false_easting (" << false_easting << ").\n\n"; - exit(1); + return; } // false_northing @@ -1820,11 +1860,11 @@ void NcCfFile::get_grid_mapping_lambert_azimuthal_equal_area(const NcVar *grid_m if(!is_bad_data(false_northing) && !is_eq(false_northing, 0.0)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "MET cannot process Lambert Azimuthal Equal Area files " << "with non-zero false_northing (" << false_northing << ").\n\n"; - exit(1); + return; } // Calculate the pin indices. The pin will be located at the grid's reference @@ -1845,6 +1885,8 @@ void NcCfFile::get_grid_mapping_lambert_azimuthal_equal_area(const NcVar *grid_m grid.set(data); if (dy_m < 0) grid.set_swap_to_north(true); + grid_ready = true; + } @@ -1863,10 +1905,10 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin grid_mapping_var, (string)"standard_parallel"); if (IS_INVALID_NC_P(std_parallel_att)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Cannot get standard_parallel attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // longitude_of_central_meridian @@ -1875,10 +1917,10 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin grid_mapping_var, (string)"longitude_of_central_meridian"); if (IS_INVALID_NC_P(central_lon_att)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Cannot get longitude_of_central_meridian attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // latitude_of_projection_origin @@ -1887,10 +1929,10 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin grid_mapping_var, (string)"latitude_of_projection_origin"); if (IS_INVALID_NC_P(proj_origin_lat_att)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Cannot get latitude_of_projection_origin attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // Look for the x/y dimensions and x/y coordinate variables @@ -1899,9 +1941,9 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin if (get_data_size(_xCoordVar) != (int) GET_NC_SIZE_P(_xDim) || get_data_size(_yCoordVar) != (int) GET_NC_SIZE_P(_yDim)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Coordinate variables don't match dimension sizes in netCDF file.\n\n"; - exit(1); + return; } // Make sure that the coordinate variables are given in meters. If we get @@ -1925,9 +1967,9 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin x_coord_units_name == "meters") x_coord_to_m_cf = 1.0; else if (x_coord_units_name == "km") x_coord_to_m_cf = 1000.0; else { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "The X coordinates must be in meters or kilometers for MET.\n\n"; - exit(1); + return; } } } @@ -1949,9 +1991,9 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin y_coord_units_name == "meters" ) y_coord_to_m_cf = 1.0; else if (y_coord_units_name == "km") y_coord_to_m_cf = 1000.0; else { - mlog << Error << "\n" << method_name << " -> " - << "The X coordinates must be in meters or kilometers for MET.\n\n"; - exit(1); + mlog << Warning << "\n" << method_name << " -> " + << "The Y coordinates must be in meters or kilometers for MET.\n\n"; + return; } } } @@ -1983,9 +2025,9 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin if (fabs(dx_m_a - dy_m_a) > DELTA_TOLERANCE) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "MET can only process Lambert Conformal files where the x-axis and y-axis deltas are the same\n\n"; - exit(1); + return; } // As a sanity check, make sure that the deltas are constant through the @@ -1996,9 +2038,9 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin double curr_delta = fabs(x_values[i] - x_values[i-1]); if (fabs(curr_delta - dx_m_a) > DELTA_TOLERANCE) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "MET can only process Lambert Conformal files where the delta along the x-axis is constant\n\n"; - exit(1); + return; } } @@ -2007,9 +2049,9 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin double curr_delta = fabs(y_values[i] - y_values[i-1]); if (fabs(curr_delta - dy_m_a) > DELTA_TOLERANCE) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "MET can only process Lambert Conformal files where the delta along the y-axis is constant\n\n"; - exit(1); + return; } } @@ -2051,6 +2093,7 @@ void NcCfFile::get_grid_mapping_lambert_conformal_conic(const NcVar *grid_mappin grid.set(data); if (dy_m < 0) grid.set_swap_to_north(true); + grid_ready = true; if(std_parallel_att) delete std_parallel_att; if(central_lon_att) delete central_lon_att; @@ -2144,32 +2187,32 @@ void NcCfFile::get_grid_mapping_latitude_longitude(const NcVar *grid_mapping_var if (_xDim == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find X dimension (degrees_east) in netCDF file.\n\n"; - exit(1); + return; } if (_yDim == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find Y dimension (degrees_north) in netCDF file.\n\n"; - exit(1); + return; } if (_xCoordVar == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find X coord variable (" << GET_NC_NAME_P(_xDim) << ") in netCDF file.\n\n"; - exit(1); + return; } if (_yCoordVar == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find Y coord variable (" << GET_NC_NAME_P(_yDim) << ") in netCDF file.\n\n"; - exit(1); + return; } long lon_counts = _xDim->getSize(); @@ -2177,9 +2220,9 @@ void NcCfFile::get_grid_mapping_latitude_longitude(const NcVar *grid_mapping_var if (get_data_size(_xCoordVar) != lon_counts || get_data_size(_yCoordVar) != lat_counts) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Coordinate variables don't match dimension sizes in netCDF file.\n\n"; - exit(1); + return; } get_grid_from_lat_lon_vars(_yCoordVar, _xCoordVar, lat_counts, lon_counts); @@ -2277,17 +2320,17 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va << "This is an ellipsoidal earth.\n\n"; } else if(!has_scale_factor && !has_standard_parallel) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "The attribute \"scale_factor_at_projection_origin\" and \"standard_parallel\" of the " << GET_NC_NAME_P(grid_mapping_var) << " variable do not exist.\n\n"; - exit(1); + return; } else if(has_scale_factor && !is_eq(proj_origin_scale_factor, 1.0)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "Unexpected attribute value of " << proj_origin_scale_factor << " for the scale_factor_at_projection_origin attribute of the " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // Look for the x/y dimensions and x/y coordinate variables @@ -2296,9 +2339,9 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va if (get_data_size(_xCoordVar) != (int) GET_NC_SIZE_P(_xDim) || get_data_size(_yCoordVar) != (int) GET_NC_SIZE_P(_yDim)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "Coordinate variables don't match dimension sizes in netCDF file.\n\n"; - exit(1); + return; } // Make sure that the coordinate variables are given in meters. If we get @@ -2322,10 +2365,10 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va x_coord_units_name == "meters") x_coord_to_m_cf = 1.0; else if ( x_coord_units_name == "km") x_coord_to_m_cf = 1000.0; else { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "The X coordinates (" << x_coord_units_name << ") must be in meters or kilometers for MET.\n\n"; - exit(1); + return; } } } @@ -2347,9 +2390,9 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va y_coord_units_name == "meters" ) y_coord_to_m_cf = 1.0; else if ( y_coord_units_name == "km") y_coord_to_m_cf = 1000.0; else { - mlog << Error << "\n" << method_name - << "The X coordinates must be in meters or kilometers for MET.\n\n"; - exit(1); + mlog << Warning << "\n" << method_name + << "The Y coordinates must be in meters or kilometers for MET.\n\n"; + return; } } } @@ -2381,9 +2424,9 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va if (fabs(dx_m_a - dy_m_a) > DELTA_TOLERANCE) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "MET can only process Polar Stereographic files where the x-axis and y-axis deltas are the same.\n\n"; - exit(1); + return; } if (is_eq(semi_major_axis, bad_data_double)) @@ -2470,6 +2513,7 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va data.dump(); grid.set(data); + grid_ready = true; //Note: do not set grid.set_swap_to_north() @@ -2611,10 +2655,10 @@ void NcCfFile::get_grid_mapping_rotated_latitude_longitude(const NcVar *grid_map grid_mapping_var, (string)"grid_north_pole_latitude"); if (IS_INVALID_NC_P(grid_np_lat_att)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Cannot get grid_north_pole_latitude attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // grid_north_pole_longitude @@ -2623,10 +2667,10 @@ void NcCfFile::get_grid_mapping_rotated_latitude_longitude(const NcVar *grid_map grid_mapping_var, (string)"grid_north_pole_longitude"); if (IS_INVALID_NC_P(grid_np_lon_att)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Cannot get grid_north_pole_longitude attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // Look for the grid_latitude and grid_longitude dimensions @@ -2703,32 +2747,32 @@ void NcCfFile::get_grid_mapping_rotated_latitude_longitude(const NcVar *grid_map if (_xDim == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find X dimension (degrees_east) in netCDF file.\n\n"; - exit(1); + return; } if (_yDim == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find Y dimension (degrees_north) in netCDF file.\n\n"; - exit(1); + return; } if (_xCoordVar == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find X coord variable (" << GET_NC_NAME_P(_xDim) << ") in netCDF file.\n\n"; - exit(1); + return; } if (_yCoordVar == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find Y coord variable (" << GET_NC_NAME_P(_yDim) << ") in netCDF file.\n\n"; - exit(1); + return; } long lon_counts = _xDim->getSize(); @@ -2736,9 +2780,9 @@ void NcCfFile::get_grid_mapping_rotated_latitude_longitude(const NcVar *grid_map if (get_data_size(_xCoordVar) != lon_counts || get_data_size(_yCoordVar) != lat_counts) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Coordinate variables don't match dimension sizes in netCDF file.\n\n"; - exit(1); + return; } // Store spacing in LatLon data structure @@ -2832,10 +2876,10 @@ void NcCfFile::get_grid_mapping_geostationary( grid_mapping_var, (string)"perspective_point_height"); if (IS_INVALID_NC_P(perspective_point_height_att)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Cannot get perspective_point_height attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // semi_major_axis @@ -2843,10 +2887,10 @@ void NcCfFile::get_grid_mapping_geostationary( grid_mapping_var, (string)"semi_major_axis"); if (IS_INVALID_NC_P(semi_major_axis_att)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Cannot get semi_major_axis attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // semi_minor_axis @@ -2854,10 +2898,10 @@ void NcCfFile::get_grid_mapping_geostationary( grid_mapping_var, (string)"semi_minor_axis"); if (IS_INVALID_NC_P(semi_minor_axis_att)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Cannot get semi_minor_axis attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // inverse_flattening @@ -2865,10 +2909,10 @@ void NcCfFile::get_grid_mapping_geostationary( grid_mapping_var, (string)"inverse_flattening"); if (IS_INVALID_NC_P(inverse_flattening_att)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Cannot get inverse_flattening attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // latitude_of_projection_origin @@ -2876,10 +2920,10 @@ void NcCfFile::get_grid_mapping_geostationary( grid_mapping_var, (string)"latitude_of_projection_origin"); if (IS_INVALID_NC_P(proj_origin_lat_att)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Cannot get latitude_of_projection_origin attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // longitude_of_projection_origin @@ -2887,10 +2931,10 @@ void NcCfFile::get_grid_mapping_geostationary( grid_mapping_var, (string)"longitude_of_projection_origin"); if (IS_INVALID_NC_P(proj_origin_lon_att)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Cannot get longitude_of_projection_origin attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // sweep_angle_axis @@ -2898,10 +2942,10 @@ void NcCfFile::get_grid_mapping_geostationary( grid_mapping_var, (string)"sweep_angle_axis"); if (IS_INVALID_NC_P(sweep_angle_axis_att)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Cannot get sweep_angle_axis attribute from " << GET_NC_NAME_P(grid_mapping_var) << " variable.\n\n"; - exit(1); + return; } // Look for the x/y dimensions and x/y coordinate variables @@ -2910,21 +2954,21 @@ void NcCfFile::get_grid_mapping_geostationary( bool do_exit = false; if (_xDim == nullptr) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Didn't find X dimension (projection_x_coordinate) in netCDF file.\n\n"; do_exit = true; } if (_yDim == nullptr) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Didn't find Y dimension (projection_y_coordinate) in netCDF file.\n\n"; do_exit = true; } if (_xCoordVar == nullptr) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Didn't find X coord variable (" << GET_NC_NAME_P(_xDim) << ") in netCDF file.\n\n"; do_exit = true; @@ -2932,7 +2976,7 @@ void NcCfFile::get_grid_mapping_geostationary( if (_yCoordVar == nullptr) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Didn't find Y coord variable (" << GET_NC_NAME_P(_yDim) << ") in netCDF file.\n\n"; do_exit = true; @@ -2941,12 +2985,12 @@ void NcCfFile::get_grid_mapping_geostationary( if (get_data_size(_xCoordVar) != (int) GET_NC_SIZE_P(_xDim) || get_data_size(_yCoordVar) != (int) GET_NC_SIZE_P(_yDim)) { - mlog << Error << "\n" << method_name + mlog << Warning << "\n" << method_name << "-> Coordinate variables don't match dimension sizes in netCDF file.\n\n"; do_exit = true; } - if (do_exit) exit(1); + if (do_exit) return; // Figure out the dx/dy and x/y pin values from the dimension variables @@ -2955,7 +2999,6 @@ void NcCfFile::get_grid_mapping_geostationary( get_nc_data(_xCoordVar, x_values); - long y_counts = GET_NC_SIZE_P(_yDim); double y_values[y_counts]; @@ -2977,7 +3020,7 @@ void NcCfFile::get_grid_mapping_geostationary( NcVar *var_y_bound = (NcVar *)nullptr; for (int j=0; j 0) { data.x_image_bounds = new double[bound_count]; data.y_image_bounds = new double[bound_count]; - if (0 != var_x_bound) get_nc_data(var_x_bound, data.x_image_bounds); - if (0 != var_y_bound) get_nc_data(var_y_bound, data.y_image_bounds); + if (nullptr != var_x_bound) get_nc_data(var_x_bound, data.x_image_bounds); + if (nullptr != var_y_bound) get_nc_data(var_y_bound, data.y_image_bounds); } double flatten = 1.0/data.inverse_flattening; @@ -3034,6 +3077,7 @@ void NcCfFile::get_grid_mapping_geostationary( // Note: Computing lat/lon was deferred because it took 1 minutes grid.set(data); + grid_ready = true; if (perspective_point_height_att) delete perspective_point_height_att; if (semi_major_axis_att) delete semi_major_axis_att; @@ -3112,19 +3156,19 @@ bool NcCfFile::get_grid_from_coordinates(const NcVar *data_var) { } if (_xCoordVar == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find X coord variable (" << x_dim_var_name << ") in netCDF file.\n\n"; if (coordinates_att) delete coordinates_att; - return true; + return false; } if (_yCoordVar == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find Y coord variable (" << y_dim_var_name << ") in netCDF file.\n\n"; if (coordinates_att) delete coordinates_att; - return true; + return false; } StringArray dimNames; @@ -3153,10 +3197,10 @@ bool NcCfFile::get_grid_from_coordinates(const NcVar *data_var) { if ((x_size != lon_counts && x_size != latlon_counts) || (y_size != lat_counts && x_size != latlon_counts)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Coordinate variables don't match dimension sizes in netCDF file.\n\n"; if (coordinates_att) delete coordinates_att; - exit(1); + return false; } if (coordinates_att) { @@ -3280,18 +3324,18 @@ bool NcCfFile::get_grid_from_dimensions() if (_xCoordVar == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find X coord variable (" << GET_NC_NAME_P(_xDim) << ") in netCDF file.\n\n"; - exit(1); + return false; } if (_yCoordVar == nullptr) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Didn't find Y coord variable (" << GET_NC_NAME_P(_yDim) << ") in netCDF file.\n\n"; - exit(1); + return false; } long lat_counts = GET_NC_SIZE_P(_yDim); @@ -3327,19 +3371,25 @@ void NcCfFile::get_grid_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, LatLonData NcCfFile::get_data_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, const long lat_counts, const long lon_counts, bool &swap_to_north) { - static const string method_name = "get_data_from_lat_lon_vars()"; + static const string method_name = "NcCfFile::get_data_from_lat_lon_vars()"; // Figure out the dlat/dlon values from the dimension variables + LatLonData data; + data.name = latlon_proj_type; + data.Nlat = (int)lat_counts; + data.Nlon = (int)lon_counts; + long x_size = get_data_size(lon_var); long y_size = get_data_size(lat_var); long latlon_counts = lon_counts*lat_counts; bool two_dim_coord = (x_size == latlon_counts) && (y_size == latlon_counts ); + if( !two_dim_coord && (x_size != lon_counts || y_size != lat_counts)) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Coordinate variables don't match dimension sizes in netCDF file.\n\n"; - exit(1); + return data; } double lat_values[lat_counts]; @@ -3365,6 +3415,8 @@ LatLonData NcCfFile::get_data_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, get_nc_data(lat_var,lat_values); get_nc_data(lon_var,lon_values); } + data.lat_ll = lat_values[0]; + data.lon_ll = rescale_lon(-lon_values[0]); // Calculate dlat and dlon assuming they are constant. @@ -3376,6 +3428,9 @@ LatLonData NcCfFile::get_data_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, << " lon[" << (lon_counts-1) << "]=" << lon_values[lon_counts-1] << " dlon=" << dlon << "\n"; + data.delta_lat = dlat; + data.delta_lon = dlon; + ConcatString point_nccf; bool skip_sanity_check = get_att_value_string(_ncFile, nc_att_met_point_nccf, point_nccf); if (!skip_sanity_check) { @@ -3406,7 +3461,7 @@ LatLonData NcCfFile::get_data_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, << i-1 << "]=" << lat_values[i-1] << " lat[" << i << "]=" << lat_values[i] << " " << fabs(curr_delta - dlat) << " > " << degree_tolerance << "\n"; - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "MET can only process Latitude/Longitude files where the latitudes are evenly spaced (dlat=" << dlat <<", delta[" << i << "]=" << curr_delta << ")\n\n"; sanity_check_failed = true; @@ -3426,7 +3481,7 @@ LatLonData NcCfFile::get_data_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, << i-1 << "]=" << lon_values[i-1] << " lon[" << i << "]=" << lon_values[i] << " " << fabs(curr_delta - dlon) << " > " << degree_tolerance << "\n"; - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "MET can only process Latitude/Longitude files where the longitudes are evenly spaced (dlon=" << dlon <<", delta[" << i << "]=" << curr_delta << ")\n\n"; sanity_check_failed = true; @@ -3435,9 +3490,9 @@ LatLonData NcCfFile::get_data_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, } if (sanity_check_failed) { - mlog << Error << "\n" << method_name << " -> " + mlog << Warning << "\n" << method_name << " -> " << "Please check the input data is the lat/lon projection\n\n"; - exit(1); + return data; } } @@ -3450,16 +3505,6 @@ LatLonData NcCfFile::get_data_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, // guaranteed anywhere that I see. But if this is not the case, then we // will probably also need to reorder the data itself. - LatLonData data; - - data.name = latlon_proj_type; - data.lat_ll = lat_values[0]; - data.lon_ll = rescale_lon(-lon_values[0]); - data.delta_lat = dlat; - data.delta_lon = dlon; - data.Nlat = lat_counts; - data.Nlon = lon_counts; - if (dlat < 0) { swap_to_north = true; data.delta_lat = -dlat; @@ -3468,6 +3513,7 @@ LatLonData NcCfFile::get_data_from_lat_lon_vars(NcVar *lat_var, NcVar *lon_var, else { swap_to_north = false; } + grid_ready = true; return data; diff --git a/src/libcode/vx_data2d_nc_cf/nc_cf_file.h b/src/libcode/vx_data2d_nc_cf/nc_cf_file.h index 5e29896459..f6f0b8b3fa 100644 --- a/src/libcode/vx_data2d_nc_cf/nc_cf_file.h +++ b/src/libcode/vx_data2d_nc_cf/nc_cf_file.h @@ -56,18 +56,18 @@ class NcCfFile { int getNx() const { - if (_xDim == 0) + if (_xDim == nullptr) return 0; - return GET_NC_SIZE_P(_xDim); + return (int)GET_NC_SIZE_P(_xDim); } int getNy() const { - if (_yDim == 0) + if (_yDim == nullptr) return 0; - return GET_NC_SIZE_P(_yDim); + return (int)GET_NC_SIZE_P(_yDim); } NcVarInfo *get_time_var_info() const { return _time_var_info; } @@ -112,6 +112,8 @@ class NcCfFile { bool getData(const char *, const LongArray &, DataPlane &, NcVarInfo *&) const; + bool update_grid(const Grid &); + Grid build_grid_from_lat_lon_vars(netCDF::NcVar *lat_var, netCDF::NcVar *lon_var, const long lat_counts, const long lon_counts); NcVarInfo* find_var_name(const char * var_name) const; @@ -122,6 +124,8 @@ class NcCfFile { static const double DELTA_TOLERANCE; netCDF::NcFile * _ncFile; // allocated + bool grid_ready; + bool has_attr_grid; // // dimensions diff --git a/src/libcode/vx_nc_util/nc_utils.hpp b/src/libcode/vx_nc_util/nc_utils.hpp index c9024a471d..7298a7e4d8 100644 --- a/src/libcode/vx_nc_util/nc_utils.hpp +++ b/src/libcode/vx_nc_util/nc_utils.hpp @@ -236,46 +236,100 @@ void apply_scale_factor_(T *data, const int cell_count, clock_t start_clock = clock(); const char *method_name = "apply_scale_factor(T) "; - if (cell_count > 0) { - int idx; - int positive_cnt = 0; - int unpacked_count = 0; - T min_value, max_value; - T raw_min_val, raw_max_val; + if (cell_count <= 0) return; + + T min_value; + T max_value; + T raw_min_val; + T raw_max_val; + int idx = 0; + int positive_cnt = 0; + int unpacked_count = 0; + + if (has_fill_value) { + // Set met_fill_value (-9999) for FillValues (missing values) + for (; idx data[idx]) raw_min_val = data[idx]; + if (raw_max_val < data[idx]) raw_max_val = data[idx]; + data[idx] = (data[idx] * scale_factor) + add_offset; + if (data[idx] > 0) positive_cnt++; + if (min_value > data[idx]) min_value = data[idx]; + if (max_value < data[idx]) max_value = data[idx]; + unpacked_count++; } + } + mlog << Debug(debug_level) << method_name << var_name + << "(data_type=" << typeid(data[0]).name() << "): unpacked data: count=" + << unpacked_count << " out of " << cell_count + << ", scale_factor=" << scale_factor<< " add_offset=" << add_offset + << ". FillValue(" << data_type << ")=" << nc_fill_value << "\n"; + mlog << Debug(debug_level) << method_name + << " data range [" << min_value << " - " << max_value + << "] raw data: [" << raw_min_val << " - " << raw_max_val + << "] Positive count: " << positive_cnt << "\n"; - raw_min_val = raw_max_val = data[idx]; - min_value = max_value = (data[idx] * scale_factor) + add_offset; - for (; idx data[idx]) raw_min_val = data[idx]; - if (raw_max_val < data[idx]) raw_max_val = data[idx]; - data[idx] = (data[idx] * scale_factor) + add_offset; - if (data[idx] > 0) positive_cnt++; - if (min_value > data[idx]) min_value = data[idx]; - if (max_value < data[idx]) max_value = data[idx]; - unpacked_count++; - } + mlog << Debug(debug_level) << method_name << " took " + << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + return; +} + +//////////////////////////////////////////////////////////////////////// + +template +void update_missing_values(T *data, const long cell_count, + const T nc_fill_value, const T met_fill_value, + const char *data_type, const char *var_name) { + int missing_count = 0; + const int debug_level = 7; + clock_t start_clock = clock(); + const char *method_name = "update_missing_values(T) "; + + if (cell_count <= 0) return; + + T max_value; + T min_value; + long idx = 0; + int positive_cnt = 0; + + // Set met_fill_value (-9999) for FillValues (missing values) + for (; idx data[idx]) min_value = data[idx]; + if (max_value < data[idx]) max_value = data[idx]; + } + } + mlog << Debug(debug_level) << method_name << var_name + << "(data_type=" << typeid(data[0]).name() << "): FillValue(" << data_type << ")=" << nc_fill_value << "\n"; + mlog << Debug(debug_level) << method_name + << " data range [" << min_value << " - " << max_value + << "] Positive count: " << positive_cnt << "\n"; + if (0 < missing_count) { + mlog << Debug(3) << method_name << var_name + << "(data_type=" << typeid(data[0]).name() << "): found " << missing_count << " FillValues out of " << cell_count << "\n"; } mlog << Debug(debug_level) << method_name << " took " << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; @@ -311,18 +365,23 @@ bool get_nc_data_(netCDF::NcVar *var, T *data, const T met_missing) { bool return_status = get_nc_data_t(var, data); if (return_status) { + T nc_missing; + const int cell_count = get_data_size(var); + bool has_missing_attr = get_var_fill_value(var, nc_missing); + if (!has_missing_attr) nc_missing = met_missing; + //scale_factor and add_offset if (has_add_offset_attr(var) || has_scale_factor_attr(var)) { - T nc_missing; - const int cell_count = get_data_size(var); double add_offset = get_var_add_offset(var); double scale_factor = get_var_scale_factor(var); - bool has_missing_attr = get_var_fill_value(var, nc_missing); - if (!has_missing_attr) nc_missing = met_missing; apply_scale_factor_(data, cell_count, add_offset, scale_factor, nc_missing, met_missing, has_missing_attr, "", GET_NC_NAME_P(var).c_str()); } + else if (has_missing_attr) { + update_missing_values(data, cell_count, nc_missing, met_missing, + "", GET_NC_NAME_P(var).c_str()); + } } return return_status; } @@ -330,7 +389,7 @@ bool get_nc_data_(netCDF::NcVar *var, T *data, const T met_missing) { //////////////////////////////////////////////////////////////////////// template -bool get_nc_data_(netCDF::NcVar *var, T *data, T bad_data, const LongArray &dims, const LongArray &curs) { +bool get_nc_data_(netCDF::NcVar *var, T *data, T met_missing, const LongArray &dims, const LongArray &curs) { bool return_status = false; const char *method_name = "get_nc_data_(T, *dims, *curs) "; @@ -371,7 +430,7 @@ bool get_nc_data_(netCDF::NcVar *var, T *data, T bad_data, const LongArray &dims } for (int idx1=0; idx1getVar(start, count, data); return_status = true; + T nc_missing; + bool has_missing_attr = get_var_fill_value(var, nc_missing); + if (!has_missing_attr) nc_missing = met_missing; + //scale_factor and add_offset if (has_add_offset_attr(var) || has_scale_factor_attr(var)) { - T nc_missing; double add_offset = get_var_add_offset(var); double scale_factor = get_var_scale_factor(var); - bool has_missing_attr = get_var_fill_value(var, nc_missing); - if (!has_missing_attr) nc_missing = bad_data; apply_scale_factor_(data, data_size, add_offset, scale_factor, - nc_missing, bad_data, has_missing_attr, + nc_missing, met_missing, has_missing_attr, "", GET_NC_NAME_P(var).c_str()); } + else if (has_missing_attr) { + update_missing_values(data, data_size, nc_missing, met_missing, + "", GET_NC_NAME_P(var).c_str()); + } } return return_status; } @@ -437,17 +501,22 @@ bool get_nc_data_(netCDF::NcVar *var, T *data, T met_missing, const long dim, co var->getVar(start, count, data); return_status = true; + T nc_missing; + bool has_missing_attr = get_var_fill_value(var, nc_missing); + if (!has_missing_attr) nc_missing = met_missing; + //scale_factor and add_offset if (has_add_offset_attr(var) || has_scale_factor_attr(var)) { - T nc_missing; double add_offset = get_var_add_offset(var); double scale_factor = get_var_scale_factor(var); - bool has_missing_attr = get_var_fill_value(var, nc_missing); - if (!has_missing_attr) nc_missing = met_missing; apply_scale_factor_(data, dim, add_offset, scale_factor, nc_missing, met_missing, has_missing_attr, "", GET_NC_NAME_P(var).c_str()); } + else if (has_missing_attr) { + update_missing_values(data, dim, nc_missing, met_missing, + "", GET_NC_NAME_P(var).c_str()); + } } return return_status; } @@ -456,7 +525,7 @@ bool get_nc_data_(netCDF::NcVar *var, T *data, T met_missing, const long dim, co // read a single data template -bool get_nc_data_(netCDF::NcVar *var, T *data, T bad_data, const LongArray &curs) { +bool get_nc_data_(netCDF::NcVar *var, T *data, T met_missing, const LongArray &curs) { bool return_status = false; //const char *method_name = "get_nc_data_(*curs) "; @@ -470,7 +539,7 @@ bool get_nc_data_(netCDF::NcVar *var, T *data, T bad_data, const LongArray &curs } // Retrieve the NetCDF value from the NetCDF variable. - return_status = get_nc_data_(var, data, bad_data, dims, curs); + return_status = get_nc_data_(var, data, met_missing, dims, curs); } return return_status; } diff --git a/src/libcode/vx_statistics/apply_mask.cc b/src/libcode/vx_statistics/apply_mask.cc index 03e31e97c7..bd12b1a25b 100644 --- a/src/libcode/vx_statistics/apply_mask.cc +++ b/src/libcode/vx_statistics/apply_mask.cc @@ -89,21 +89,8 @@ Grid parse_grid_string(const char *grid_str) { Grid grid; StringArray sa; - // Parse as a white-space separated string - sa.parse_wsss(grid_str); - - // Search for a named grid - if(sa.n() == 1 && find_grid_by_name(sa[0].c_str(), grid)) { - mlog << Debug(3) << "Use the grid named \"" - << grid_str << "\".\n"; - } - // Parse grid definition - else if(sa.n() > 1 && parse_grid_def(sa, grid)) { - mlog << Debug(3) << "Use the grid defined by string \"" - << grid_str << "\".\n"; - } - // Extract the grid from a gridded data file - else { + if (!build_grid_by_grid_string(grid_str, grid, "parse_grid_strin", false)) { + // Extract the grid from a gridded data file mlog << Debug(3) << "Use the grid defined by file \"" << grid_str << "\".\n"; @@ -228,24 +215,8 @@ void parse_grid_mask(const ConcatString &mask_grid_str, Grid &grid) { // Check for empty input string if(mask_grid_str.empty()) return; - // Parse mask_grid_str as a white-space separated string - StringArray sa; - sa.parse_wsss(mask_grid_str); - - // Named grid - if(sa.n() == 1 && find_grid_by_name(mask_grid_str.c_str(), grid)) { - mlog << Debug(3) - << "Use the grid named \"" << mask_grid_str << "\".\n"; - } - // Grid specification string - else if(sa.n() > 1 && parse_grid_def(sa, grid)) { - mlog << Debug(3) - << "Use the grid defined by string \"" << mask_grid_str - << "\".\n"; - } - // Extract the grid from a gridded data file - else { - + if (!build_grid_by_grid_string(mask_grid_str, grid, "parse_grid_mask", false)) { + // Extract the grid from a gridded data file mlog << Debug(3) << "Use the grid defined by file \"" << mask_grid_str << "\".\n"; diff --git a/src/tools/other/gen_vx_mask/gen_vx_mask.cc b/src/tools/other/gen_vx_mask/gen_vx_mask.cc index 7612f0d7d2..98b15a212a 100644 --- a/src/tools/other/gen_vx_mask/gen_vx_mask.cc +++ b/src/tools/other/gen_vx_mask/gen_vx_mask.cc @@ -177,24 +177,8 @@ void process_command_line(int argc, char **argv) { void process_input_grid(DataPlane &dp) { - // Parse the input grid as a white-space separated string - StringArray sa; - sa.parse_wsss(input_gridname); - - // Search for a named grid - if(sa.n() == 1 && find_grid_by_name(sa[0].c_str(), grid)) { - mlog << Debug(3) - << "Use input grid named \"" << input_gridname << "\".\n"; - } - // Parse grid definition - else if(sa.n() > 1 && parse_grid_def(sa, grid)) { - mlog << Debug(3) - << "Use input grid defined by string \"" << input_gridname - << "\".\n"; - } - // Extract the grid from a gridded data file - else { - + if (!build_grid_by_grid_string(input_gridname, grid, "process_input_grid", false)) { + // Extract the grid from a gridded data file mlog << Debug(3) << "Use input grid defined by file \"" << input_gridname << "\".\n"; @@ -284,22 +268,7 @@ void process_mask_file(DataPlane &dp) { // For the grid mask type, support named grids and grid // specification strings if(mask_type == MaskType::Grid) { - - // Parse the mask file as a white-space separated string - StringArray sa; - sa.parse_wsss(mask_filename); - - // Search for a named grid - if(sa.n() == 1 && find_grid_by_name(sa[0].c_str(), grid_mask)) { - mlog << Debug(3) - << "Use mask grid named \"" << mask_filename << "\".\n"; - } - // Parse grid definition - else if(sa.n() > 1 && parse_grid_def(sa, grid_mask)) { - mlog << Debug(3) - << "Use mask grid defined by string \"" << mask_filename - << "\".\n"; - } + build_grid_by_grid_string(mask_filename, grid_mask, "process_mask_file", false); } // Parse as a gridded data file if not already set diff --git a/src/tools/other/mode_time_domain/3d_conv.cc b/src/tools/other/mode_time_domain/3d_conv.cc index d204ee61b0..1ecbe0b68a 100644 --- a/src/tools/other/mode_time_domain/3d_conv.cc +++ b/src/tools/other/mode_time_domain/3d_conv.cc @@ -154,8 +154,6 @@ struct DataHandle { { - int j; - _out.put('\n'); _out << " DataHandle:\n"; @@ -168,13 +166,13 @@ struct DataHandle { _out << " plane_loaded = ["; - for (j=0; j 0) { - - // Parse as a white-space separated string - sa.parse_wsss(plot_grid_string); - - // Search for a named grid - if(sa.n() == 1 && find_grid_by_name(sa[0].c_str(), grid)) { - mlog << Debug(3) << "Use the grid named \"" - << plot_grid_string << "\".\n"; - } - // Parse grid definition - else if(sa.n() > 1 && parse_grid_def(sa, grid)) { - mlog << Debug(3) << "Use the grid defined by string \"" - << plot_grid_string << "\".\n"; - } - // Extract the grid from a gridded data file - else { + if (!build_grid_by_grid_string(plot_grid_string, grid, + "PlotPointObsConfInfo::process_config -> ", false)) { + // Extract the grid from a gridded data file mlog << Debug(3) << "Use the grid defined by file \"" << plot_grid_string << "\".\n"; diff --git a/src/tools/other/point2grid/point2grid.cc b/src/tools/other/point2grid/point2grid.cc index de7d284fc1..73210e4b27 100644 --- a/src/tools/other/point2grid/point2grid.cc +++ b/src/tools/other/point2grid/point2grid.cc @@ -30,7 +30,7 @@ #include #include -#include +#include #include "main.h" #include "vx_log.h" @@ -424,6 +424,7 @@ static void process_data_file() { // Get the obs type before opening NetCDF obs_type = get_obs_type(nc_in); goes_data = (obs_type == TYPE_GOES || obs_type == TYPE_GOES_ADP); + if (obs_type == TYPE_UNKNOWN && ftype == FileType_NcCF) obs_type = TYPE_NCCF; if (obs_type == TYPE_NCCF) setenv(nc_att_met_point_nccf, "yes", 1); // Read the input data file @@ -594,14 +595,33 @@ bool get_nc_data_string_array(NcFile *nc, const char *var_name, static int get_obs_type(NcFile *nc) { int obs_type = TYPE_UNKNOWN; + MetConfig config; ConcatString att_val_scene_id; ConcatString att_val_project; ConcatString input_type; static const char *method_name = "get_obs_type() -> "; - bool has_project = get_global_att(nc, (string)"project", att_val_project); - bool has_scene_id = get_global_att(nc, (string)"scene_id", att_val_scene_id); - if( has_scene_id && has_project && att_val_project == "GOES" ) { + bool has_attr_grid = false; + auto vinfo = VarInfoFactory::new_var_info(FileType_NcCF); + for(int i=0; iclear(); + // Populate the VarInfo object using the config string + config.read_string(FieldSA[i].c_str()); + vinfo->set_dict(config); + if (vinfo->grid_attr().is_set()) { + has_attr_grid = true; + break; + } + } + if (vinfo) { delete vinfo; vinfo = (VarInfo *) nullptr; } + + if (has_attr_grid) { + obs_type = TYPE_NCCF; + input_type = "OBS_NCCF"; + } + else if (get_global_att(nc, (string)"scene_id", att_val_scene_id) + && get_global_att(nc, (string)"project", att_val_project) + && att_val_project == "GOES" ) { obs_type = TYPE_GOES; input_type = "GOES"; if (!adp_filename.empty()) { @@ -693,12 +713,15 @@ std::set prepare_qoes_qc_array() { void process_point_met_data(MetPointData *met_point_obs, MetConfig &config, VarInfo *vinfo, const Grid &to_grid) { - int idx, hdr_idx; + int hdr_idx; int var_idx_or_gc; ConcatString vname; - DataPlane fr_dp, to_dp; - DataPlane cnt_dp, mask_dp; - DataPlane prob_dp, prob_mask_dp; + DataPlane fr_dp; + DataPlane to_dp; + DataPlane cnt_dp; + DataPlane mask_dp; + DataPlane prob_dp; + DataPlane prob_mask_dp; bool has_prob_thresh = !prob_cat_thresh.check(bad_data_double); @@ -807,7 +830,7 @@ void process_point_met_data(MetPointData *met_point_obs, MetConfig &config, VarI } else { bool not_found_grib_code = true; - for (idx=0; idxobs_ids[idx]) { not_found_grib_code = false; break; @@ -824,7 +847,7 @@ void process_point_met_data(MetPointData *met_point_obs, MetConfig &config, VarI if (exit_by_field_name_error) { ConcatString log_msg; if (use_var_id) { - for (idx=0; idxobs_ids[idx])) { grib_codes.add(obs_data->obs_ids[idx]); if (0 < idx) log_msg << ", "; @@ -855,7 +878,9 @@ void process_point_met_data(MetPointData *met_point_obs, MetConfig &config, VarI // Check the time range. Apply the time window bool valid_time_from_config = true; - unixtime valid_beg_ut, valid_end_ut, obs_time; + unixtime valid_beg_ut; + unixtime valid_end_ut; + unixtime obs_time; valid_time_array.clear(); valid_time = vinfo->valid(); @@ -865,7 +890,7 @@ void process_point_met_data(MetPointData *met_point_obs, MetConfig &config, VarI valid_beg_ut = valid_end_ut = valid_time; if (!is_bad_data(conf_info.beg_ds)) valid_beg_ut += conf_info.beg_ds; if (!is_bad_data(conf_info.end_ds)) valid_end_ut += conf_info.end_ds; - for(idx=0; idx valid_time) valid_time = obs_time; } @@ -915,7 +940,7 @@ void process_point_met_data(MetPointData *met_point_obs, MetConfig &config, VarI int filtered_by_time = 0; int filtered_by_msg_type = 0; int filtered_by_qc = 0; - for (idx=0; idx < nobs; idx++) { + for (int idx=0; idx < nobs; idx++) { if (var_idx_or_gc == obs_data->obs_ids[idx]) { var_count2++; hdr_idx = obs_data->obs_hids[idx]; @@ -1104,7 +1129,7 @@ void process_point_met_data(MetPointData *met_point_obs, MetConfig &config, VarI log_msg << ", by msg_type: " << filtered_by_msg_type; if (0 < filtered_by_msg_type) { log_msg << " ["; - for(idx=0; idx 0) log_msg << ","; log_msg << conf_info.message_type[idx]; } @@ -1113,7 +1138,7 @@ void process_point_met_data(MetPointData *met_point_obs, MetConfig &config, VarI log_msg << ", by QC: " << filtered_by_qc; if (0 < filtered_by_qc) { log_msg << " ["; - for(idx=0; idx 0) log_msg << ","; log_msg << qc_flags[idx]; } @@ -1179,7 +1204,7 @@ static void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, nc_point_obs.close(); mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; } @@ -1216,7 +1241,7 @@ static void process_point_python(const string python_command, MetConfig &config, met_point_file.close(); mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; return; } @@ -1234,10 +1259,10 @@ static void process_point_nccf_file(NcFile *nc_in, MetConfig &config, bool *skip_times = nullptr; double *valid_times = nullptr; int filtered_by_time = 0; + int time_from_size = 1; clock_t start_clock = clock(); bool opt_all_attrs = false; Grid fr_grid = fr_mtddf->grid(); - int from_size = fr_grid.nx() * fr_grid.ny(); static const char *method_name = "process_point_nccf_file() -> "; NcVar var_lat; @@ -1247,10 +1272,10 @@ static void process_point_nccf_file(NcFile *nc_in, MetConfig &config, ConcatString lon_vname = conf_info.get_var_name(conf_key_lon_vname); if (lat_vname != conf_key_lat_vname && lon_vname != conf_key_lon_vname) { - if (lat_vname != conf_key_lat_vname && has_var(nc_in, lat_vname.c_str())) { + if (has_var(nc_in, lat_vname.c_str())) { var_lat = get_nc_var(nc_in, lat_vname.c_str()); } - if (lon_vname != conf_key_lon_vname && has_var(nc_in, lon_vname.c_str())) { + if (has_var(nc_in, lon_vname.c_str())) { var_lon = get_nc_var(nc_in, lon_vname.c_str()); } if (IS_INVALID_NC(var_lat)) { @@ -1270,7 +1295,7 @@ static void process_point_nccf_file(NcFile *nc_in, MetConfig &config, // Find lat/lon variables from the coordinates attribue if (0 < FieldSA.n() && !user_defined_latlon) { ConcatString coordinates_value; - VarInfoNcCF var_info = VarInfoNcCF(*(VarInfoNcCF *)vinfo); + auto var_info = VarInfoNcCF(*(VarInfoNcCF *)vinfo); // Initialize var_info.clear(); // Populate the VarInfo object using the config string @@ -1322,8 +1347,9 @@ static void process_point_nccf_file(NcFile *nc_in, MetConfig &config, if( IS_VALID_NC(time_var) ) { if( 1 < get_dim_count(&time_var) ) { double max_time = bad_data_double; - skip_times = new bool[from_size]; - valid_times = new double[from_size]; + time_from_size = get_data_size(&time_var); + skip_times = new bool[time_from_size]; + valid_times = new double[time_from_size]; if (get_nc_data(&time_var, valid_times)) { int sec_per_unit = 0; bool no_leap_year = false; @@ -1334,7 +1360,7 @@ static void process_point_nccf_file(NcFile *nc_in, MetConfig &config, if (!is_bad_data(conf_info.end_ds)) valid_end_ut += conf_info.end_ds; ref_ut = get_reference_unixtime(&time_var, sec_per_unit, no_leap_year); } - for (int i=0; i 0 ) { tmp_time = add_to_unixtime(ref_ut, sec_per_unit, valid_times[i], no_leap_year); @@ -1405,7 +1431,10 @@ static void process_point_nccf_file(NcFile *nc_in, MetConfig &config, // List range of data values if(mlog.verbosity_level() >= 2) { - double fr_dmin, fr_dmax, to_dmin, to_dmax; + double fr_dmin; + double fr_dmax; + double to_dmin; + double to_dmax; fr_dp.data_range(fr_dmin, fr_dmax); to_dp.data_range(to_dmin, to_dmax); mlog << Debug(2) << "Range of data (" << FieldSA[i] << ")\n" @@ -1428,7 +1457,8 @@ static void process_point_nccf_file(NcFile *nc_in, MetConfig &config, bool has_prob_thresh = !prob_cat_thresh.check(bad_data_double); if (has_prob_thresh || do_gaussian_filter) { - DataPlane prob_dp, prob_mask_dp; + DataPlane prob_dp; + DataPlane prob_mask_dp; ConcatString vname_prob = vname; vname_prob << "_prob_" << prob_cat_thresh.get_abbr_str(); int nx = to_dp.nx(); @@ -1477,12 +1507,12 @@ static void process_point_nccf_file(NcFile *nc_in, MetConfig &config, cellMapping = (IntArray *) nullptr; if( 0 < filtered_by_time ) { mlog << Debug(2) << method_name << "Filtered by time: " << filtered_by_time - << " out of " << from_size + << " out of " << time_from_size << " [" << unix_to_yyyymmdd_hhmmss(valid_beg_ut) << " to " << unix_to_yyyymmdd_hhmmss(valid_end_ut) << "]\n"; } mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; return; } @@ -1506,9 +1536,6 @@ static void regrid_nc_variable(NcFile *nc_in, Met2dDataFile *fr_mtddf, exit(1); } - int from_lat_cnt = fr_grid.ny(); - int from_lon_cnt = fr_grid.nx(); - int from_data_size = from_lat_cnt * from_lon_cnt; if(!fr_mtddf->data_plane(*vinfo, fr_dp)) { mlog << Error << "\n" << method_name << "Trouble reading data \"" @@ -1516,106 +1543,112 @@ static void regrid_nc_variable(NcFile *nc_in, Met2dDataFile *fr_mtddf, << InputFilename << "\"\n\n"; exit(1); } - else { - bool is_to_north = !fr_grid.get_swap_to_north(); - auto from_data = new float[from_data_size]; - for (int xIdx=0; xIdx= 4) { - if (from_min_value > data_value) from_min_value = data_value; - if (from_max_value < data_value) from_max_value = data_value; - } + int from_index; + int no_map_cnt = 0; + int missing_cnt = 0; + int non_missing_cnt = 0; + double data_value; + IntArray cellArray; + NumArray dataArray; + double from_min_value = 10e10; + double from_max_value = -10e10; + int to_lat_cnt = to_grid.ny(); + int to_lon_cnt = to_grid.nx(); + + missing_cnt = non_missing_cnt = 0; + to_dp.set_constant(bad_data_double); + + for (int xIdx=0; xIdx= 9) { - double to_lat, to_lon; - to_grid.xy_to_latlon(xIdx,yIdx, to_lat, to_lon); - to_lon *= -1; - if (1 == data_cnt) - mlog << Debug(9) << method_name - << "value: " << to_value << " to (" << to_lon << ", " << to_lat - << ") from offset " << from_index << ".\n"; - else - mlog << Debug(9) << method_name - << "value: " << to_value - << ", max: " << dataArray.max() - << ", min: " << dataArray.min() - << ", mean: " << dataArray.sum()/data_cnt - << " from " << data_cnt << " (out of " << cellArray.n() - << ") data values to (" << to_lon << ", " << to_lat << ").\n"; - } + dataArray.add(data_value); + non_missing_cnt++; + if(mlog.verbosity_level() >= 4) { + if (from_min_value > data_value) from_min_value = data_value; + if (from_max_value < data_value) from_max_value = data_value; } } - else { - no_map_cnt++; + + if (0 < dataArray.n()) { + double to_value; + int data_cnt = dataArray.n(); + if (1 == data_cnt) to_value = dataArray[0]; + else if (RGInfo.method == InterpMthd::Min) to_value = dataArray.min(); + else if (RGInfo.method == InterpMthd::Max) to_value = dataArray.max(); + else if (RGInfo.method == InterpMthd::Median) { + dataArray.sort_array(); + to_value = dataArray[data_cnt/2]; + if (0 == data_cnt % 2) + to_value = (to_value + dataArray[(data_cnt/2)+1])/2; + } + else to_value = dataArray.sum() / data_cnt; // UW_Mean + + to_dp.set(to_value, xIdx, yIdx); + to_cell_cnt++; + if(mlog.verbosity_level() >= 9) { + double to_lat; + double to_lon; + to_grid.xy_to_latlon(xIdx,yIdx, to_lat, to_lon); + to_lon *= -1; + if (1 == data_cnt) + mlog << Debug(9) << method_name + << "value: " << to_value << " to (" << to_lon << ", " << to_lat + << ") from offset " << from_index << ".\n"; + else + mlog << Debug(9) << method_name + << "value: " << to_value + << ", max: " << dataArray.max() + << ", min: " << dataArray.min() + << ", mean: " << dataArray.sum()/data_cnt + << " from " << data_cnt << " (out of " << cellArray.n() + << ") data values to (" << to_lon << ", " << to_lat << ").\n"; + } } } + else { + no_map_cnt++; + } } + } - delete [] from_data; + delete [] from_data; - mlog << Debug(4) << method_name << "[Count] data cells: " << to_cell_cnt - << ", missing: " << missing_cnt << ", non_missing: " << non_missing_cnt - << ", non mapped cells: " << no_map_cnt - << " out of " << (to_lat_cnt*to_lon_cnt) - << "\n\tRange: data: [" << from_min_value << " - " << from_max_value - << "]\n"; - } + mlog << Debug(4) << method_name << "[Count] data cells: " << to_cell_cnt + << ", missing: " << missing_cnt << ", non_missing: " << non_missing_cnt + << ", non mapped cells: " << no_map_cnt + << " out of " << (to_lat_cnt*to_lon_cnt) + << "\n\tRange: data: [" << from_min_value << " - " << from_max_value + << "]\n"; if (to_cell_cnt == 0) { mlog << Warning << "\n" << method_name @@ -1623,7 +1656,7 @@ static void regrid_nc_variable(NcFile *nc_in, Met2dDataFile *fr_mtddf, } mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; } //////////////////////////////////////////////////////////////////////// @@ -1818,7 +1851,10 @@ static void process_goes_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, // List range of data values if(mlog.verbosity_level() >= 2) { - double fr_dmin, fr_dmax, to_dmin, to_dmax; + double fr_dmin; + double fr_dmax; + double to_dmin; + double to_dmax; fr_dp.data_range(fr_dmin, fr_dmax); to_dp.data_range(to_dmin, to_dmax); mlog << Debug(2) << "Range of data (" << FieldSA[i] << ")\n" @@ -1848,7 +1884,8 @@ static void process_goes_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, bool has_prob_thresh = !prob_cat_thresh.check(bad_data_double); if (has_prob_thresh || do_gaussian_filter) { - DataPlane prob_dp, prob_mask_dp; + DataPlane prob_dp; + DataPlane prob_mask_dp; ConcatString vname_prob = vname; vname_prob << "_prob_" << prob_cat_thresh.get_abbr_str(); int nx = to_dp.nx(); @@ -1901,7 +1938,7 @@ static void process_goes_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, clear_cell_mapping(cellMapping); delete [] cellMapping; cellMapping = (IntArray *) nullptr; mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; return; } @@ -2018,7 +2055,8 @@ static bool get_grid_mapping(const Grid &to_grid, IntArray *cellMapping, return status; } - double x, y; + double x; + double y; DataPlane to_dp; int to_lat_count = to_grid.ny(); int to_lon_count = to_grid.nx(); @@ -2050,7 +2088,7 @@ static bool get_grid_mapping(const Grid &to_grid, IntArray *cellMapping, << ((obs_count > 0) ? 1.0*count_in_grid/obs_count*100 : 0) << "%)\n"; } mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; return status; } @@ -2061,15 +2099,17 @@ static void get_grid_mapping_latlon( DataPlane from_dp, DataPlane to_dp, Grid to_grid, IntArray *cellMapping, float *latitudes, float *longitudes, int from_lat_count, int from_lon_count, bool *skip_times, bool to_north) { - double x, y; - double to_ll_lat, to_ll_lon; + double x; + double y; + double to_ll_lat; + double to_ll_lon; int count_in_grid = 0; clock_t start_clock = clock(); int to_lat_count = to_grid.ny(); int to_lon_count = to_grid.nx(); int to_size = to_lat_count * to_lon_count; int data_size = from_lat_count * from_lon_count; - static const char *method_name = "get_grid_mapping(lats, lons) -> "; + static const char *method_name = "get_grid_mapping_latlon(lats, lons) -> "; auto to_cell_counts = new int[to_size]; auto mapping_indices = new int[data_size]; @@ -2088,7 +2128,7 @@ static void get_grid_mapping_latlon( float lat = latitudes[coord_offset]; float lon = longitudes[coord_offset]; if( lat < MISSING_LATLON || lon < MISSING_LATLON ) continue; - to_grid.latlon_to_xy(lat, -1.0*lon, x, y); + to_grid.latlon_to_xy(lat, -1.0*rescale_lon(lon), x, y); int idx_x = nint(x); int idx_y = nint(y); if (0 <= idx_x && idx_x < to_lon_count && 0 <= idx_y && idx_y < to_lat_count) { @@ -2097,16 +2137,17 @@ static void get_grid_mapping_latlon( to_cell_counts[to_offset] += 1; count_in_grid++; if(mlog.verbosity_level() >= 15) { - double to_lat, to_lon; + double to_lat; + double to_lon; to_grid.xy_to_latlon(idx_x, idx_y, to_lat, to_lon); - mlog << Debug(15) << method_name << " [" << xIdx << "," << yIdx << "] to " << coord_offset - << " (" << lon << ", " << lat << ") to (" << (to_lon*-1) << ", " << to_lat << ")\n"; + mlog << Debug(15) << method_name << " index: [" << xIdx << "," << yIdx << "] to " << coord_offset + << " (" << lon << ", " << lat << ") to (" << rescale_lon(-1*to_lon) << ", " << to_lat << ")\n"; } } } } mlog << Debug(LEVEL_FOR_PERFORMANCE+2) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds for mapping cells\n"; + << get_exe_duration(start_clock) << " seconds for mapping cells\n"; // Count the mapping cells for each to_cell and prepare IntArray int max_count = 0; @@ -2129,7 +2170,7 @@ static void get_grid_mapping_latlon( } } mlog << Debug(LEVEL_FOR_PERFORMANCE+1) << method_name << "took " - << (clock()-tmp_clock)/double(CLOCKS_PER_SEC) + << get_exe_duration(tmp_clock) << " seconds for extending IntArray (max_cells=" << max_count << ")\n"; // Build cell mapping @@ -2150,7 +2191,7 @@ static void get_grid_mapping_latlon( mlog << Debug(3) << method_name << "within grid: " << count_in_grid << " out of " << data_size << " (" << 1.0*count_in_grid/data_size*100 << "%)\n"; mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; } //////////////////////////////////////////////////////////////////////// @@ -2158,7 +2199,8 @@ static void get_grid_mapping_latlon( static bool get_grid_mapping(const Grid &fr_grid, const Grid &to_grid, IntArray *cellMapping, NcVar var_lat, NcVar var_lon, bool *skip_times) { bool status = false; - DataPlane from_dp, to_dp; + DataPlane from_dp; + DataPlane to_dp; ConcatString cur_coord_name; clock_t start_clock = clock(); static const char *method_name = "get_grid_mapping(var_lat, var_lon) -> "; @@ -2167,7 +2209,14 @@ static bool get_grid_mapping(const Grid &fr_grid, const Grid &to_grid, IntArray int to_lon_count = to_grid.nx(); int from_lat_count = fr_grid.ny(); int from_lon_count = fr_grid.nx(); - + if (0 == from_lat_count) { + int dim_offset = get_dim_count(&var_lat) - 1; + if (dim_offset < 0) dim_offset = 0; + from_lat_count = get_dim_size(&var_lat, dim_offset); + dim_offset = get_dim_count(&var_lon) - 2; + if (dim_offset < 0) dim_offset = 0; + from_lon_count = get_dim_size(&var_lon, dim_offset); + } // Override the from nx & ny from NetCDF if exists int data_size = from_lat_count * from_lon_count; mlog << Debug(4) << method_name << "data_size (ny*nx): " << data_size @@ -2212,7 +2261,7 @@ static bool get_grid_mapping(const Grid &fr_grid, const Grid &to_grid, IntArray if( longitudes ) delete [] longitudes; } // if data_size > 0 mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; return status; } @@ -2272,7 +2321,8 @@ static ConcatString get_goes_grid_input(MetConfig config, const Grid fr_grid) { static void get_grid_mapping(const Grid &fr_grid, const Grid &to_grid, IntArray *cellMapping, const ConcatString &geostationary_file) { static const char *method_name = "get_grid_mapping() -> "; - DataPlane from_dp, to_dp; + DataPlane from_dp; + DataPlane to_dp; ConcatString cur_coord_name; clock_t start_clock = clock(); @@ -2431,7 +2481,7 @@ static void get_grid_mapping(const Grid &fr_grid, const Grid &to_grid, IntArray if(coord_nc_in) delete coord_nc_in; mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; } //////////////////////////////////////////////////////////////////////// @@ -2831,7 +2881,7 @@ static void regrid_goes_variable(NcFile *nc_in, const VarInfo *vinfo, } mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; } //////////////////////////////////////////////////////////////////////// @@ -2900,7 +2950,7 @@ static void save_geostationary_data(const ConcatString geostationary_file, delete nc_file; nc_file = nullptr; mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " - << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; + << get_exe_duration(start_clock) << " seconds\n"; } ////////////////////////////////////////////////////////////////////////