From 2a2c1c4775e828121e84dfa050bca25f1aef839f Mon Sep 17 00:00:00 2001 From: johnhg Date: Tue, 26 Jan 2021 16:37:03 -0700 Subject: [PATCH] Update develop-ref after #1639 (#1642) * Getting rid of compiler warnings in PB2NC by replacing several instances of the NULL pointer with the nul character (\0) instead. * Fix typo in config_options.rst. * Feature 1408 var_name_for_grib_code (#1617) * #1408 Added get_var_id * #1408 Check variable name in the configuration to use the variable name instewad of grib code * #1408 Added point2grid_ascii2nc_surfrad_DW_PSP_by_name * Feature 1580 2d time (#1616) * #1580 Added get_grid_from_lat_lon_vars * #1580 Added get_grid_from_lat_lon_vars and support 2D time variable * #1580 Support int type variable without scale_factor and add_offset attributes * #1580 Support 2D time variable. Implemented filtering by valid_time * #1580 Bug fix: read time with dimension 0 * #1580 Support time variable with no dimension * #1580 Initial release * #1580 Added point2grid_2D_time * #1580 Check project attribute for GOES * #1580 Changed NULL to 0 to avoid co,pilation warning * #1580 Added point2grid_2D_time * #1580 Added "point2grid configuration file" section * #1580 Changed to_grid for point2grid_NCCF_UK & point2grid_2D_time Co-authored-by: Howard Soh Co-authored-by: John Halley Gotway * feature 1580 nccf (#1619) * #1580 Correct the precision at _apply_scale_factor * #1580 Added unit test plot_data_plane_NCCF_time * #1580 Changed argument type to double at _apply_scale_factor(double) * Bugfix 1618 develop pb2nc (#1623) Co-authored-by: Howard Soh * Feature 1624 OBS_COMMAND (#1625) * Per #1627, add grid_data.regrid config option for PlotPointObs and update the tool to do the requested regridding. Still need to update the docs. * Per #1627, update docs about grid_data.regrid config option for PlotPointObs. * Per #1627, add another call to plot_point_obs to exercise the new regrid functionality. * Feature 1624 obs_command second try (#1629) * Per #1624, define OBS_COMMAND. * Per #1624, unset the test-specific environment variables after completing the run. * Per #1624, after PR #1625 merged these changes into develop, they caused 2 unexpected diffs in the NB output. These were caused by enviornment variables being unset after each test. Updating unit_netcdf.xml and unit_point2grid.xml to define more test-specific environment variables to reproduce previous NB output. * Organizing NB climatology and point2grid output files into the appopriate directories rather than having them at the top-level directory. * Update pull_request_template.md * Update the point2grid unit tests to write their temp files to the point2grid subdirectory instead of the top-level test output directory. * Update appendixC.rst Split the definition of H_RATE and POD * Feature 1626 tc_gen (#1633) * Per #1448, many changes for TC-Gen. Replace the oper_genesis dictionary with the oper_technique string. Add genesis_init_diff config entry. Update config_constants.h accordingly and the tc_gen_conf_info.h/.cc to parse the updated config entries. * Per #1448, large overhaul of the tc_gen matching logic. This work is not yet complete. Still need to compute categorical MISSES but the current version does compile. * Per #1448, add GenesisInfoArray::has_storm_id() function and remove the unused set_dland() function. * Per #1448, more updates. Define the best genesis events while parsing the best tracks. We need to know the best genesis events in order to count up the forecast misses. * Per #1448, lots more changes for tc_gen. Create a PairDataGenesis class to store genesis pairs. This will be needed to write a matched pair line type. * Per #1448, minor tweaks to log messages. * Per #1448, update PairDataGenesis class to store the BEST track Storm ID since the forecast genesis do not have meaningful Storm ID's. * Per #1448, in GenesisInfoArray::add(), do NOT store multiple genesis events for the same storm, but do print a useful Debug(3) log message about it. * Per #1448, update PairDataGenesis::has_case() logic to check the storm id and initialization time but NOT require an exact forecast hour match. * Per #1448, update the tc_gen log messages to more concisely and consistently report the storm id. * Per #1448, update the PairDataGenesis logic a bit to have all the misses and hits in chronological order. * Per #1448, add genesis_init_diff entry. * Per #1448, set the default genesis_init_diff entry to 48 hours since that's what Dan H used in his examples. * Per #1448, work on comments and log messages. * Per #1448, reimplement TrackInfoArray as a vector instead of managing the memory myself. This makes the implmentation of TrackInfoArray::erase_storm_id() very easy. Replace n_tracks() function with n() in several places. * Per #1448, add valid_freq and basin_file config entries. Also rename load_dland.h/.cc to load_tc_data.h/.cc and add code to read the basin file. * Per #1448, add GenesisInfoArray::erase_storm_id(). * Per #1448, update tc_gen code to handle new config options. * Per #1448, had my units wrong. Was processing seconds when I thought it was hours! * Per #1448, making test TC-Gen config file consistent with the default. * Per #1448, also track the obs valid times. * Per #1448, switch from tech1/tech2 to dev/ops methods. Update log messages and add lots of details to the tc_gen documentation. * Per #1430, in tc_gen enable dev_method_flag, ops_method_flag, ci_alpha, and output_flag to be specified separately for each filter. Also add nc_pairs_flag and genesis_track_points_window config options. Add config constants entries for these options and update tc_gen to handle all of these changes. * Per #1430, consolidate the parse_grid_mask() code a bit to avoid redundancy.: * Per #1430, just cleaning up some messy comments. * Per #1430, adding hooks for writing NetCDF output file. * Per #1430, update DataPlane::set_size() function to take a 3rd argument to specify how the DataPlane should be initialized. * Per #1430, update the nc_pairs_flag options and update the code to parse them. * Per #1430, update the TrackInfo class to track and report the min/max warm core information. * Per #1430, current state of development. Still a work in progress. I'm getting runtime segfaults when testing and I still need to NOT overcount the BEST track hits. * Per #1430, committing changes described by https://github.com/dtcenter/MET/issues/1430#issuecomment-755755898 * Per #1430, forgot to rename genesis_match_window to genesis_hit_window as it is in the code. * Per #1430, chaning GenesisInfo to just inherit directly from TrackInfo. Frankly, I should have thought of this a LONG time ago. * Per #1430, change the default desc setting from NA to ALL and add the best_unique_flag option. * Per #1430, simplify the logic now that GenesisInfo is derived from TrackInfo. Also support the best_unique_flag config option. * Per #1430, instead of storing 12 individual DataPlane objects, store them in a map to make writing their output more convenient. * Per #1430, updating documentation and comments. * Per #1430, more doc updates. * Per #1430, update unit test to only write NetCDF counts for the AL_BASIN and not the other filters. * Per #1430, fix parsing logic for nc_pairs_flag = TRUE. * Per #1430, fix bug. Check the VxOpt.NcInfo before calling write_nc(), not the top-level one. * Per #1430, the docker build of tc_gen failed. * Per #1430, working on DockerHub compilation. * Per #1430, getting DockerHub build working. * One more try. * Per #1597, add hooks for new GENMPR stat line type. * Per #1597, add config file option and column definitions for the GENMPR line type. * Per #1597, finish writing the GENMPR line type. * Per #1597, change the default output grid from a global 5 degree to global 1 degree grid. * Per #1597, change GENMPR output columns to GEN_TDIFF and INIT_TDIFF since they're reported in HHMMSS format instead of seconds. Also, tweak the config file for the tc-gen unit test. * Per #1597, have to add GENMPR header columns for Stat-Analysis and test scripts to handle it. * Per #1597, update Stat-Analysis to handle the GENMPR line type. * Per #1597, user's guide updates for the GENMPR and NetCDF output file. * Per #1597, add AGEN_INIT and AGEN_FHR columns. * Per #1597, add AGEN_INIT and AGEN_FHR columns. * Per #1597, remove the AGEN_TIME and BGEN_TIME columns from the GENMPR line type and instead write the genesis times to the FCST_VALID_BEG/END and OBS_VALID_BEG/END header columns. * Remove some unused output column name definitions. There are a remnant from very early versions of MET which included the CTP, CFP, and COP line types. * Per #1597, update config file options to use dev_hit_radius, dev_hit_window, and opt_hit_tdiff. Also update log message to switch from 'lead' to 'forecast hour'. * Per #1626, add met_regrid_nearest() utility function since I'm calling it twice. * Per #1626, update the basin_global_tenth_degree.nc basin definition file to include basin name abbreviations. * Per #1626, update load_tc_data.h/.cc to also read the basin abbreviations from the NetCDF basin file. * Per #1626, add TC-Gen config file options for init_inc, init_exc, and basin_mask. Updated the library and application code, and updated the user's guide. * Fixing Fortify warnings for 'Poor Style: Variable Never Used' in 6 files. * Fix Fortify warnings for 'Uninitialized variable' in tc_gen.cc and point2grid.cc. * Fix Fortify warnings for 'Poor Style: Redundant Initialization' in plot_point_obs.cc and point2grid.cc. * Feature 1346 valid time attr (#1634) * #1346 get_att_value_unixtime supports yyyymmdd_hhmmss, too * #1346 Check valid_time & init_time attributes, too * #1346 Check valid_time & init_time attributes, too Co-authored-by: Howard Soh * Feature 1473 python errors (#1615) * Added sample script to read ascii data and create an xarray. * Disabled use_xarray exit for testing. * Get attrs from DataArray if using xarray. * Removed some comments. * Revised error messages for use with both numpy and xarray. * Removing commented out code. Co-authored-by: David Fillmore Co-authored-by: johnhg * Feature 1630 zero obs (#1637) * Per #1630, update ascii2nc to change zero observations from an error (which returns bad status) to a warning message. * Per #1630, update point2grid to read an empty input file and write fields of 0's or bad data to the output. Change previous error message to warning. Also, update LOTS of warning and error log messages to make them consistent. * Per #1630, need to initialize the dataplanes before the loop (for when there are no obs) and within each loop iteration (for when there are multiple fields to process). * Bugfix 1638 develop climo cdf (#1639) * Per #1638, correct the order of arguments in the call to the normal_cdf() utility function. * Per #1638, update the logic in derive_climo_prob(). For CDP thresholds, the constant climo probability should be based on the inequality type where less-than-types match the threshold percentile value while greater-than-types are 1.0 minus the threshold percentile. * Per #1638, update normal_cdf() to initialize the output CDF field using the climo mean field instead of the observation data field. This makes the timestamps consistent for the climo mean, stdev, and cdf variables in the Grid-Stat NetCDF matched pairs output file. * Update tc_gen.cc Co-authored-by: hsoh-u Co-authored-by: Howard Soh Co-authored-by: John Halley Gotway Co-authored-by: j-opatz <59586397+j-opatz@users.noreply.github.com> Co-authored-by: David Fillmore Co-authored-by: David Fillmore --- met/src/basic/vx_util/data_plane_util.cc | 2 +- met/src/libcode/vx_data2d_nc_met/met_file.cc | 14 +- .../dataplane_from_numpy_array.cc | 8 +- met/src/libcode/vx_nc_util/load_tc_data.cc | 1 - met/src/libcode/vx_nc_util/nc_var_info.cc | 45 ++-- met/src/libcode/vx_statistics/pair_base.cc | 35 ++- met/src/libcode/vx_tc_util/genesis_info.cc | 1 - met/src/tools/core/grid_stat/grid_stat.cc | 8 +- met/src/tools/other/ascii2nc/ascii2nc.cc | 1 + met/src/tools/other/ascii2nc/file_handler.cc | 8 +- met/src/tools/other/grid_diag/grid_diag.cc | 185 +++++++------ .../other/plot_point_obs/plot_point_obs.cc | 7 +- met/src/tools/other/point2grid/point2grid.cc | 248 ++++++++++-------- met/src/tools/tc_utils/tc_gen/tc_gen.cc | 117 ++++++++- .../tools/tc_utils/tc_gen/tc_gen_conf_info.cc | 4 +- 15 files changed, 423 insertions(+), 261 deletions(-) diff --git a/met/src/basic/vx_util/data_plane_util.cc b/met/src/basic/vx_util/data_plane_util.cc index 766764fdbb..7174a36b84 100644 --- a/met/src/basic/vx_util/data_plane_util.cc +++ b/met/src/basic/vx_util/data_plane_util.cc @@ -606,7 +606,7 @@ DataPlane subtract(const DataPlane &dp1, const DataPlane &dp2) { DataPlane normal_cdf(const DataPlane &dp, const DataPlane &mn, const DataPlane &sd) { - DataPlane cdf = dp; + DataPlane cdf = mn; double v; // Check grid dimensions diff --git a/met/src/libcode/vx_data2d_nc_met/met_file.cc b/met/src/libcode/vx_data2d_nc_met/met_file.cc index 8efadbf5ee..6f922b80e6 100644 --- a/met/src/libcode/vx_data2d_nc_met/met_file.cc +++ b/met/src/libcode/vx_data2d_nc_met/met_file.cc @@ -35,8 +35,10 @@ using namespace std; static const char x_dim_name [] = "lon"; static const char y_dim_name [] = "lat"; -static const string valid_time_att_name = "valid_time_ut"; -static const string init_time_att_name = "init_time_ut"; +static const string valid_time_att_name = "valid_time"; +static const string init_time_att_name = "init_time"; +static const string valid_time_ut_att_name = "valid_time_ut"; +static const string init_time_ut_att_name = "init_time_ut"; static const string accum_time_att_name = "accum_time_sec"; static const string name_att_name = "name"; @@ -295,8 +297,12 @@ for (j=0; j " - << "numpy array is not 2-dimensional! ... " + mlog << Error << "\ndataplane_from_data_array() -> " + << "data array is not 2-dimensional! ... " << "(dim = " << (np.n_dims()) << ")\n\n"; exit ( 1 ); @@ -170,8 +170,8 @@ else if ( dtype == ">f8" ) load_numpy (np.buffer(), Nx, Ny, bi else { - mlog << Error << "\ndataplane_from_numpy_array() -> " - << "unsupported numpy data type \"" << dtype << "\"\n\n"; + mlog << Error << "\ndataplane_from_data_array() -> " + << "unsupported data type \"" << dtype << "\"\n\n"; exit ( 1 ); diff --git a/met/src/libcode/vx_nc_util/load_tc_data.cc b/met/src/libcode/vx_nc_util/load_tc_data.cc index 39710f912e..240ac92c01 100644 --- a/met/src/libcode/vx_nc_util/load_tc_data.cc +++ b/met/src/libcode/vx_nc_util/load_tc_data.cc @@ -86,7 +86,6 @@ void load_tc_basin(const ConcatString &basin_file, Grid &grid, ConcatString file_name; LongArray dim; NcVarInfo *vi; - int i; // Get the path for the distance to land file file_name = replace_path(basin_file); diff --git a/met/src/libcode/vx_nc_util/nc_var_info.cc b/met/src/libcode/vx_nc_util/nc_var_info.cc index 28609fba9e..416e2e2c4a 100644 --- a/met/src/libcode/vx_nc_util/nc_var_info.cc +++ b/met/src/libcode/vx_nc_util/nc_var_info.cc @@ -28,6 +28,7 @@ using namespace std; #include "vx_cal.h" unixtime get_att_value_unixtime(const NcAtt *att) { + ConcatString s; unixtime time_value = -1; switch ( GET_NC_TYPE_ID_P(att) ) { case NC_INT64: @@ -36,13 +37,18 @@ unixtime get_att_value_unixtime(const NcAtt *att) { break; case NC_CHAR: - ConcatString s; get_att_value_chars(att, s); - time_value = string_to_unixtime(s.c_str()); + // 20120410_120000 VS. 1333929600 + if (0 > s.find('_') && 0 > s.find('-')) + time_value = string_to_unixtime(s.c_str()); + else + time_value = yyyymmdd_hhmmss_to_unix(s.c_str()); break; - //default: - // break; + default: + mlog << Warning << "get_att_value_unixtime() The attribute type (" + << GET_NC_TYPE_NAME_P(att) << ") is not supported\n"; + break; } // switch return time_value; } @@ -326,12 +332,11 @@ bool get_att_int(const NcVarInfo &info, const ConcatString att_name, int &att_va { - NcVarAtt *att; - bool found = false; att_value = bad_data_int; - att = get_nc_att(info.var, att_name, false); - if (!IS_INVALID_NC_P(att)) { + NcVarAtt *att = get_nc_att(info.var, att_name, false); + bool found = IS_VALID_NC_P(att); + if (found) { att_value = get_att_value_int(att); // Check for the correct type @@ -344,7 +349,6 @@ bool get_att_int(const NcVarInfo &info, const ConcatString att_name, int &att_va exit ( 1 ); } - found = true; } if (att) delete att; @@ -364,24 +368,11 @@ bool get_att_unixtime(const NcVarInfo &info, const ConcatString att_name, unixti { - NcVarAtt *att; - bool found = false; - att_value = (unixtime) bad_data_int; - - - att = get_nc_att(info.var, att_name, false); - if (!IS_INVALID_NC_P(att)) { - found = true; - } - - if ( !found ) { - if (att) delete att; - return ( false ); - } - - // Check the type - att_value = get_att_value_unixtime(att); + + NcVarAtt *att = get_nc_att(info.var, att_name, false); + bool found = IS_VALID_NC_P(att); + if( found ) att_value = get_att_value_unixtime(att); if (att) delete att; @@ -389,7 +380,7 @@ bool get_att_unixtime(const NcVarInfo &info, const ConcatString att_name, unixti // done // - return ( true ); + return ( found ); } diff --git a/met/src/libcode/vx_statistics/pair_base.cc b/met/src/libcode/vx_statistics/pair_base.cc index c4e5a85272..bfcebaad0c 100644 --- a/met/src/libcode/vx_statistics/pair_base.cc +++ b/met/src/libcode/vx_statistics/pair_base.cc @@ -1026,9 +1026,40 @@ NumArray derive_climo_prob(const NumArray &mn_na, const NumArray &sd_na, n_mn = mn_na.n_valid(); n_sd = sd_na.n_valid(); - // For CDP threshold types, the climo_prob is constant. + // For CDP threshold types, the climo probability is constant if(othresh.get_ptype() == perc_thresh_climo_dist) { - climo_prob.add_const(othresh.get_pvalue()/100.0, n_mn); + + // Climo probability varies based on the threshold type + switch(othresh.get_type()) { + + case thresh_lt: + case thresh_le: + prob = othresh.get_pvalue()/100.0; + break; + + case thresh_eq: + prob = 0.0; + break; + + case thresh_ne: + prob = 1.0; + break; + + case thresh_gt: + case thresh_ge: + prob = 1.0 - othresh.get_pvalue()/100.0; + break; + + default: + mlog << Error << "\n\nderive_climo_prob() -> " + << "climatological threshold \"" << othresh.get_str() + << "\" cannot be converted to a probability!\n\n"; + exit(1); + break; + } + + // Add constant climo probability value + climo_prob.add_const(prob, n_mn); } // If both mean and standard deviation were provided, use them to // derive normal climatological probabilities for the current event diff --git a/met/src/libcode/vx_tc_util/genesis_info.cc b/met/src/libcode/vx_tc_util/genesis_info.cc index 218d33770b..8cc5ea36d9 100644 --- a/met/src/libcode/vx_tc_util/genesis_info.cc +++ b/met/src/libcode/vx_tc_util/genesis_info.cc @@ -442,7 +442,6 @@ const GenesisInfo & GenesisInfoArray::operator[](int n) const { //////////////////////////////////////////////////////////////////////// bool GenesisInfoArray::add(const GenesisInfo &gi) { - int i; // Skip true duplicates if(has(gi)) { diff --git a/met/src/tools/core/grid_stat/grid_stat.cc b/met/src/tools/core/grid_stat/grid_stat.cc index 46e579f765..5a90a6464f 100644 --- a/met/src/tools/core/grid_stat/grid_stat.cc +++ b/met/src/tools/core/grid_stat/grid_stat.cc @@ -1070,7 +1070,7 @@ void process_scores() { conf_info.vx_opt[i].interp_info.field); } if(conf_info.vx_opt[i].nc_info.do_climo && !cmn_dp.is_empty() && !csd_dp.is_empty()) { - write_nc((string)"CLIMO_CDF", normal_cdf(cmn_dp, csd_dp, obs_dp), + write_nc((string)"CLIMO_CDF", normal_cdf(obs_dp, cmn_dp, csd_dp), i, mthd, pnts, conf_info.vx_opt[i].interp_info.field); } @@ -1803,11 +1803,11 @@ void get_mask_points(const MaskPlane &mask_mp, apply_mask(*obs_ptr, mask_mp, pd.o_na); pd.n_obs = pd.o_na.n(); - if(cmn_ptr) apply_mask(*cmn_ptr, mask_mp, pd.cmn_na); + if(cmn_ptr) apply_mask(*cmn_ptr, mask_mp, pd.cmn_na); else pd.cmn_na.add_const(bad_data_double, pd.n_obs); - if(csd_ptr) apply_mask(*csd_ptr, mask_mp, pd.csd_na); + if(csd_ptr) apply_mask(*csd_ptr, mask_mp, pd.csd_na); else pd.csd_na.add_const(bad_data_double, pd.n_obs); - if(wgt_ptr) apply_mask(*wgt_ptr, mask_mp, pd.wgt_na); + if(wgt_ptr) apply_mask(*wgt_ptr, mask_mp, pd.wgt_na); else pd.wgt_na.add_const(default_grid_weight, pd.n_obs); if(cmn_ptr && csd_ptr) pd.add_climo_cdf(); diff --git a/met/src/tools/other/ascii2nc/ascii2nc.cc b/met/src/tools/other/ascii2nc/ascii2nc.cc index 08dcb9c1a1..4ced5397b1 100644 --- a/met/src/tools/other/ascii2nc/ascii2nc.cc +++ b/met/src/tools/other/ascii2nc/ascii2nc.cc @@ -42,6 +42,7 @@ // 014 07-23-18 Halley Gotway Support masks defined by gen_vx_mask. // 015 03-20-19 Fillmore Add aeronetv2 and aeronetv3 options. // 016 01-30-20 Bullock Add python option. +// 017 01-25-21 Halley Gotway MET #1630 Handle zero obs. // //////////////////////////////////////////////////////////////////////// diff --git a/met/src/tools/other/ascii2nc/file_handler.cc b/met/src/tools/other/ascii2nc/file_handler.cc index ceefd6529f..28ef53e8cc 100644 --- a/met/src/tools/other/ascii2nc/file_handler.cc +++ b/met/src/tools/other/ascii2nc/file_handler.cc @@ -131,14 +131,10 @@ bool FileHandler::writeNetcdfFile(const string &nc_filename) _countHeaders(); - // Check for no data + // Print a warning for no data if (_nhdr == 0) { - mlog << Error << "\nZero observations retained.\n" - << "Cannot create NetCDF Observation file: " - << nc_filename << "\n\n"; - - return false; + mlog << Warning << "\nZero observations retained!\n\n"; } mlog << Debug(2) << "Processing observations for " << _nhdr diff --git a/met/src/tools/other/grid_diag/grid_diag.cc b/met/src/tools/other/grid_diag/grid_diag.cc index 5a958847e1..2fe311e0b2 100644 --- a/met/src/tools/other/grid_diag/grid_diag.cc +++ b/met/src/tools/other/grid_diag/grid_diag.cc @@ -402,9 +402,9 @@ void setup_histograms(void) { bin_max.clear(); bin_mid.clear(); for(int k=0; kmagic_str()] = bin_min; @@ -414,9 +414,9 @@ void setup_histograms(void) { // Initialize histograms mlog << Debug(2) - << "Initializing " << data_info->magic_str() - << " histogram with " << n_bins << " bins from " - << min << " to " << max << ".\n"; + << "Initializing " << data_info->magic_str() + << " histogram with " << n_bins << " bins from " + << min << " to " << max << ".\n"; histograms[data_info->magic_str()] = vector(); init_pdf(n_bins, histograms[data_info->magic_str()]); } // for i_var @@ -426,103 +426,102 @@ void setup_histograms(void) { void setup_joint_histograms(void) { - for(int i_var=0; i_varn_bins(); + VarInfo *data_info = conf_info.data_info[i_var]; + int n_bins = data_info->n_bins(); - for(int j_var=i_var+1; j_varn_bins(); + VarInfo *joint_info = conf_info.data_info[j_var]; + int n_joint_bins = joint_info->n_bins(); - ConcatString joint_str = data_info->magic_str(); - joint_str.add("_"); - joint_str.add(joint_info->magic_str()); - mlog << Debug(2) - << "Initializing " << joint_str << " joint histogram with " - << n_bins << " x " << n_joint_bins << " bins.\n"; - joint_histograms[joint_str] = vector(); + ConcatString joint_str = data_info->magic_str(); + joint_str.add("_"); + joint_str.add(joint_info->magic_str()); + mlog << Debug(2) + << "Initializing " << joint_str << " joint histogram with " + << n_bins << " x " << n_joint_bins << " bins.\n"; + joint_histograms[joint_str] = vector(); - init_joint_pdf(n_bins, n_joint_bins, - joint_histograms[joint_str]); - } // end for j_var - } // end for i_var + init_joint_pdf(n_bins, n_joint_bins, + joint_histograms[joint_str]); + } // end for j_var + } // end for i_var } //////////////////////////////////////////////////////////////////////// void setup_nc_file(void) { - int n; - ConcatString cs; - - // Create NetCDF file - nc_out = open_ncfile(out_file.c_str(), true); - - if(IS_INVALID_NC_P(nc_out)) { - mlog << Error << "\nsetup_nc_file() -> " - << "trouble opening output NetCDF file " - << out_file << "\n\n"; - exit(1); - } - - // Add global attributes - write_netcdf_global(nc_out, out_file.c_str(), program_name, - NULL, NULL, conf_info.desc.c_str()); - add_att(nc_out, "mask_grid", (conf_info.mask_grid_name.nonempty() ? - (string)conf_info.mask_grid_name : - na_str)); - add_att(nc_out, "mask_poly", (conf_info.mask_poly_name.nonempty() ? - (string)conf_info.mask_poly_name : - na_str)); - - // Add time range information to the global attributes - add_att(nc_out, "init_beg", (string)unix_to_yyyymmdd_hhmmss(init_beg)); - add_att(nc_out, "init_end", (string)unix_to_yyyymmdd_hhmmss(init_end)); - add_att(nc_out, "valid_beg", (string)unix_to_yyyymmdd_hhmmss(valid_beg)); - add_att(nc_out, "valid_end", (string)unix_to_yyyymmdd_hhmmss(valid_end)); - add_att(nc_out, "lead_beg", (string)sec_to_hhmmss(lead_beg)); - add_att(nc_out, "lead_end", (string)sec_to_hhmmss(lead_end)); - - // Write the grid size, mask size, and series length - write_nc_var_int("grid_size", "number of grid points", grid.nxy()); - write_nc_var_int("mask_size", "number of mask points", conf_info.mask_area.count()); - write_nc_var_int("n_series", "length of series", n_series); - - // Compression level - int deflate_level = compress_level; - if(deflate_level < 0) deflate_level = conf_info.conf.nc_compression(); - - for(int i_var=0; i_var < conf_info.get_n_data(); i_var++) { - - VarInfo *data_info = conf_info.data_info[i_var]; - - // Set variable NetCDF name - ConcatString var_name = data_info->name_attr(); - var_name.add("_"); - var_name.add(data_info->level_attr()); - - // Define histogram dimensions - NcDim var_dim = add_dim(nc_out, var_name, - (long) data_info->n_bins()); - data_var_dims.push_back(var_dim); - - // Define histogram bins - ConcatString var_min_name = var_name; - ConcatString var_max_name = var_name; - ConcatString var_mid_name = var_name; - var_min_name.add("_min"); - var_max_name.add("_max"); - var_mid_name.add("_mid"); - NcVar var_min = add_var(nc_out, var_min_name, ncFloat, var_dim, - deflate_level); - NcVar var_max = add_var(nc_out, var_max_name, ncFloat, var_dim, - deflate_level); - NcVar var_mid = add_var(nc_out, var_mid_name, ncFloat, var_dim, - deflate_level); - - // Add variable attributes - cs << cs_erase << "Minimum value of " << var_name << " bin"; + ConcatString cs; + + // Create NetCDF file + nc_out = open_ncfile(out_file.c_str(), true); + + if(IS_INVALID_NC_P(nc_out)) { + mlog << Error << "\nsetup_nc_file() -> " + << "trouble opening output NetCDF file " + << out_file << "\n\n"; + exit(1); + } + + // Add global attributes + write_netcdf_global(nc_out, out_file.c_str(), program_name, + NULL, NULL, conf_info.desc.c_str()); + add_att(nc_out, "mask_grid", (conf_info.mask_grid_name.nonempty() ? + (string)conf_info.mask_grid_name : + na_str)); + add_att(nc_out, "mask_poly", (conf_info.mask_poly_name.nonempty() ? + (string)conf_info.mask_poly_name : + na_str)); + + // Add time range information to the global attributes + add_att(nc_out, "init_beg", (string)unix_to_yyyymmdd_hhmmss(init_beg)); + add_att(nc_out, "init_end", (string)unix_to_yyyymmdd_hhmmss(init_end)); + add_att(nc_out, "valid_beg", (string)unix_to_yyyymmdd_hhmmss(valid_beg)); + add_att(nc_out, "valid_end", (string)unix_to_yyyymmdd_hhmmss(valid_end)); + add_att(nc_out, "lead_beg", (string)sec_to_hhmmss(lead_beg)); + add_att(nc_out, "lead_end", (string)sec_to_hhmmss(lead_end)); + + // Write the grid size, mask size, and series length + write_nc_var_int("grid_size", "number of grid points", grid.nxy()); + write_nc_var_int("mask_size", "number of mask points", conf_info.mask_area.count()); + write_nc_var_int("n_series", "length of series", n_series); + + // Compression level + int deflate_level = compress_level; + if(deflate_level < 0) deflate_level = conf_info.conf.nc_compression(); + + for(int i_var=0; i_var < conf_info.get_n_data(); i_var++) { + + VarInfo *data_info = conf_info.data_info[i_var]; + + // Set variable NetCDF name + ConcatString var_name = data_info->name_attr(); + var_name.add("_"); + var_name.add(data_info->level_attr()); + + // Define histogram dimensions + NcDim var_dim = add_dim(nc_out, var_name, + (long) data_info->n_bins()); + data_var_dims.push_back(var_dim); + + // Define histogram bins + ConcatString var_min_name = var_name; + ConcatString var_max_name = var_name; + ConcatString var_mid_name = var_name; + var_min_name.add("_min"); + var_max_name.add("_max"); + var_mid_name.add("_mid"); + NcVar var_min = add_var(nc_out, var_min_name, ncFloat, var_dim, + deflate_level); + NcVar var_max = add_var(nc_out, var_max_name, ncFloat, var_dim, + deflate_level); + NcVar var_mid = add_var(nc_out, var_mid_name, ncFloat, var_dim, + deflate_level); + + // Add variable attributes + cs << cs_erase << "Minimum value of " << var_name << " bin"; add_var_att_local(&var_min, "long_name", cs); add_var_att_local(&var_min, "units", data_info->units_attr()); diff --git a/met/src/tools/other/plot_point_obs/plot_point_obs.cc b/met/src/tools/other/plot_point_obs/plot_point_obs.cc index 00b401d3f1..768f7393a1 100644 --- a/met/src/tools/other/plot_point_obs/plot_point_obs.cc +++ b/met/src/tools/other/plot_point_obs/plot_point_obs.cc @@ -147,7 +147,7 @@ int main(int argc, char *argv[]) { //////////////////////////////////////////////////////////////////////// void process_point_obs(const char *point_obs_filename) { - int i, h, v; + int h, v; int obs_hid_block[DEF_NC_BUFFER_SIZE]; int obs_vid_block[DEF_NC_BUFFER_SIZE]; int obs_qty_block[DEF_NC_BUFFER_SIZE]; @@ -230,10 +230,6 @@ void process_point_obs(const char *point_obs_filename) { // Get the corresponding header: // message type, staton_id, valid_time, and lat/lon/elv NcHeaderData header_data = get_nc_hdr_data(obsVars); - int strl_len = header_data.strl_len; - int typ_len = header_data.typ_len; - int sid_len = header_data.sid_len; - int vld_len = header_data.vld_len; bool use_hdr_arr = !IS_INVALID_NC(obsVars.hdr_arr_var); bool use_obs_arr = !IS_INVALID_NC(obsVars.obs_arr_var); @@ -381,7 +377,6 @@ void process_point_obs(const char *point_obs_filename) { int typ_idx, sid_idx, vld_idx; for(int i_offset=0; i_offsetname(); if (var_names.has(vname)) { mlog << Error << "\n" << method_name - << "Please add -name argument to avoid the output name conflicts on handling the variable \"" + << "Please add -name argument to avoid the output name " + << "conflicts on handling the variable \"" << vname << "\" multiple times.\n\n"; usage(); } @@ -333,7 +331,7 @@ void process_data_file() { GrdFileType ftype; ConcatString run_cs; NcFile *nc_in = (NcFile *)0; - static const char *method_name = "process_data_file() "; + static const char *method_name = "process_data_file() -> "; // Initialize configuration object MetConfig config; @@ -350,7 +348,7 @@ void process_data_file() { mlog << Debug(1) << "Reading data file: " << InputFilename << "\n"; nc_in = open_ncfile(InputFilename.c_str()); - //Get the obs type before opening NetCDF + // Get the obs type before opening NetCDF int obs_type = get_obs_type(nc_in); bool goes_data = (obs_type == TYPE_GOES || obs_type == TYPE_GOES_ADP); @@ -362,7 +360,7 @@ void process_data_file() { fr_mtddf = m_factory.new_met_2d_data_file(InputFilename.c_str(), ftype); if(!fr_mtddf) { - mlog << Error << "\n" << method_name << " -> " + mlog << Error << "\n" << method_name << "\"" << InputFilename << "\" not a valid data file\n\n"; exit(1); } @@ -376,7 +374,7 @@ void process_data_file() { vinfo = v_factory.new_var_info(ftype); if(!vinfo) { - mlog << Error << "\n" << method_name << " -> " + mlog << Error << "\n" << method_name << "unable to determine file type of \"" << InputFilename << "\"\n\n"; exit(1); @@ -417,7 +415,8 @@ void process_data_file() { } else { mlog << Error << "\n" << method_name - << "Please check the input file. Only supports GOES, MET point obs. and CF complaint NetCDF with time/lat/lon variables.\n\n"; + << "Please check the input file. Only supports GOES, MET point obs, " + << "and CF complaint NetCDF with time/lat/lon variables.\n\n"; exit(1); } @@ -529,7 +528,8 @@ int get_obs_type(NcFile *nc) { obs_type = TYPE_GOES_ADP; input_type = "GOES_ADP"; if (!file_exists(AdpFilename.c_str())) { - mlog << Error << method_name << "ADP input \"" << AdpFilename << "\" does not exist!\n"; + mlog << Error << "\n" << method_name + << "ADP input \"" << AdpFilename << "\" does not exist!\n\n"; exit(1); } } @@ -569,8 +569,9 @@ void prepare_message_types(const StringArray hdr_types) { msg_list_str << message_type_str; } if (0 == message_type_list.n()) { - mlog << Error << method_name << " No matching message types [" - << msg_list_str << "].\n"; + mlog << Error << "\n" << method_name + << "No matching message types [" + << msg_list_str << "].\n\n"; exit(1); } } @@ -637,7 +638,6 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, int *obs_hids = new int[nobs]; float *hdr_lats = new float[nhdr]; float *hdr_lons = new float[nhdr]; - //float *hdr_elvs[nhdr] = new float[nhdr]; float *obs_lvls = new float[nobs]; float *obs_hgts = new float[nobs]; float *obs_vals = new float[nobs]; @@ -683,7 +683,6 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, success_to_read = get_nc_data_float_array(nc_in, nc_var_hdr_lat, hdr_lats); if (success_to_read) success_to_read = get_nc_data_float_array(nc_in, nc_var_hdr_lon, hdr_lons); - //if (!get_nc_data_float_array(nc_in, nc_var_hdr_elv, hdr_elvs)) exit(1); if (success_to_read) success_to_read = get_nc_data_float_array(nc_in, nc_var_obs_lvl, obs_lvls); if (success_to_read) @@ -700,18 +699,22 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, bool has_qc_flags = (qc_flags.n() > 0); IntArray qc_idx_array = prepare_qc_array(qc_flags, qc_tables); - // Initialize configuration object - + // Initialize size and values of output fields nx = to_grid.nx(); ny = to_grid.ny(); to_dp.set_size(nx, ny); + to_dp.set_constant(bad_data_double); cnt_dp.set_size(nx, ny); + cnt_dp.set_constant(0); mask_dp.set_size(nx, ny); + mask_dp.set_constant(0); if (has_prob_thresh || do_gaussian_filter) { prob_dp.set_size(nx, ny); + prob_dp.set_constant(0); prob_mask_dp.set_size(nx, ny); + prob_mask_dp.set_constant(0); } - + // Loop through the requested fields int obs_count_zero_to, obs_count_non_zero_to; int obs_count_zero_from, obs_count_non_zero_from; @@ -797,7 +800,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, mlog << Error << "\n" << method_name << error_msg << "Try setting the \"name\" in the \"-field\" command line option to one of the available names:\n" - << "\t" << log_msg <<"\n\n"; + << "\t" << log_msg << "\n\n"; exit(1); } @@ -832,7 +835,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, if (obs_time > valid_time) valid_time = obs_time; } } - mlog << Debug(3) << method_name << " valid_time from " + mlog << Debug(3) << method_name << "valid_time from " << (valid_time_from_config ? "config" : "input data") << ": " << unix_to_yyyymmdd_hhmmss(valid_time) << "\n"; @@ -855,7 +858,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, ? conf_info.get_var_name(vinfo->name()) : conf_info.get_var_name(VarNameSA[i]); mlog << Debug(4) << method_name - << " var: " << vname << ", index: " << var_idx_or_gc << ".\n"; + << "var: " << vname << ", index: " << var_idx_or_gc << ".\n"; var_count = var_count2 = to_count = 0; filtered_by_time = filtered_by_msg_type = filtered_by_qc = 0; @@ -907,6 +910,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, float from_min_value = 10e10; float from_max_value = -10e10; + // Initialize counter and output fields to_count = 0; to_dp.set_constant(bad_data_double); cnt_dp.set_constant(0); @@ -984,7 +988,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, << ", mean: " << dataArray.sum()/data_count << " from " << data_count << " data values.\n"; } - mlog << Debug(8) << method_name << " data at " << x_idx << "," << y_idx + mlog << Debug(8) << method_name << "data at " << x_idx << "," << y_idx << ", value: " << to_value << "\n"; } } @@ -1042,7 +1046,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, } vinfo->set_long_name(var_long_name.c_str()); - mlog << Debug(7) << method_name << " obs_count_zero_to: " << obs_count_zero_to + mlog << Debug(7) << method_name << "obs_count_zero_to: " << obs_count_zero_to << ", obs_count_non_zero_to: " << obs_count_non_zero_to << "\n"; ConcatString log_msg; @@ -1072,22 +1076,20 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, int filtered_count = filtered_by_msg_type + filtered_by_qc + requested_valid_time; if (0 == var_count) { if (0 == filtered_count) { - mlog << Error << method_name << " No valid data for the variable [" - << vinfo->name() << "]\n"; + mlog << Warning << "\n" << method_name + << "No valid data for the variable [" + << vinfo->name() << "]\n\n"; } else { - mlog << Error << method_name << " No valid data after filtering.\n\t" - << log_msg << ".\n"; + mlog << Warning << "\n" << method_name + << "No valid data after filtering.\n\t" + << log_msg << ".\n\n"; } - exit(1); } else { - mlog << Debug(2) << method_name << " var_count=" << var_count + mlog << Debug(2) << method_name << "var_count=" << var_count << ", grid: " << to_count << " out of " << (nx * ny) << " " - << (0 < filtered_count ? log_msg.c_str() : " ") - //<< ", Obs_value] zero: " << obs_count_zero_from - //<< ", non_zero: " < 0 ) { if (!is_eq(bad_data_int, conf_info.beg_ds)) valid_beg_ut += conf_info.beg_ds; if (!is_eq(bad_data_int, conf_info.end_ds)) valid_end_ut += conf_info.end_ds; @@ -1294,7 +1297,7 @@ void regrid_nc_variable(NcFile *nc_in, Met2dDataFile *fr_mtddf, int to_cell_cnt = 0; clock_t start_clock = clock(); Grid fr_grid = fr_mtddf->grid(); - static const char *method_name = "regrid_nc_variable() --> "; + static const char *method_name = "regrid_nc_variable() -> "; NcVar var_data = get_nc_var(nc_in, vinfo->name().c_str()); if (IS_INVALID_NC(var_data)) { @@ -1308,8 +1311,10 @@ void regrid_nc_variable(NcFile *nc_in, Met2dDataFile *fr_mtddf, 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 \"" - << vinfo->name() << "\" from file \"" << InputFilename << "\"\n\n"; + mlog << Error << "\n" << method_name + << "Trouble reading data \"" + << vinfo->name() << "\" from file \"" + << InputFilename << "\"\n\n"; exit(1); } else { @@ -1396,7 +1401,7 @@ void regrid_nc_variable(NcFile *nc_in, Met2dDataFile *fr_mtddf, delete [] from_data; - mlog << Debug(4) << method_name << " Count] data cells: " << to_cell_cnt + 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) @@ -1565,7 +1570,7 @@ void process_goes_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, bool opt_all_attrs = false; clock_t start_clock = clock(); NcFile *nc_adp = (NcFile *)0; - static const char *method_name = "process_goes_file() "; + static const char *method_name = "process_goes_file() -> "; ConcatString tmp_dir = config.get_tmp_dir(); ConcatString geostationary_file(tmp_dir); @@ -1576,7 +1581,8 @@ void process_goes_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, if (0 < AdpFilename.length() && file_exists(AdpFilename.c_str())) { nc_adp = open_ncfile(AdpFilename.c_str()); if (IS_INVALID_NC_P(nc_adp)) { - mlog << Error << method_name << "Can't open the ADP input \"" << AdpFilename << "\"\n"; + mlog << Error << "\n" << method_name + << "Can't open the ADP input \"" << AdpFilename << "\"\n\n"; exit(1); } else if (is_time_mismatch(nc_in, nc_adp)) { @@ -1711,11 +1717,11 @@ void check_lat_lon(int data_size, float *latitudes, float *longitudes) { float max_lat=-90.0; float min_lon=360.0; float max_lon=-360.0; - static const char *method_name = "check_lat_lon() "; + static const char *method_name = "check_lat_lon() -> "; for (int idx=0; idx MISSING_LATLON && longitudes[idx] > MISSING_LATLON) { - mlog << Debug(7) << method_name << " index: " << idx << " lat: " + mlog << Debug(7) << method_name << "index: " << idx << " lat: " << latitudes[idx] << ", lon: " << longitudes[idx] << "\n"; cnt_printed++; } @@ -1735,15 +1741,15 @@ void check_lat_lon(int data_size, float *latitudes, float *longitudes) { mlog << Debug(7) << method_name << "Count: missing -> lat: " << cnt_missing_lat << ", lon: " << cnt_missing_lon << " invalid -> lat: " << cnt_bad_lat << ", lon: " << cnt_bad_lon << "\n" - << method_name << " LAT: " << min_lat << " to " << max_lat << "\n" - << method_name << " LONG: " << min_lon << " to " << max_lon << "\n"; + << method_name << "LAT: " << min_lat << " to " << max_lat << "\n" + << method_name << "LONG: " << min_lon << " to " << max_lon << "\n"; } //////////////////////////////////////////////////////////////////////// static unixtime compute_unixtime(NcVar *time_var, unixtime var_value) { unixtime obs_time = bad_data_int; - static const char *method_name = "compute_unixtime() "; + static const char *method_name = "compute_unixtime() -> "; if (IS_VALID_NC_P(time_var)) { int sec_per_unit; @@ -1764,11 +1770,12 @@ static bool get_grid_mapping(Grid to_grid, IntArray *cellMapping, const float *hdr_lats, const float *hdr_lons) { bool status = false; clock_t start_clock = clock(); - static const char *method_name = "get_grid_mapping(MET_obs) "; + static const char *method_name = "get_grid_mapping(MET_obs) -> "; int obs_count = obs_index_array.n(); if (0 == obs_count) { - mlog << Warning << method_name << " no valid point observation data\n"; + mlog << Warning << "\n" << method_name + << "no valid point observation data!\n\n"; return status; } @@ -1800,10 +1807,11 @@ static bool get_grid_mapping(Grid to_grid, IntArray *cellMapping, } if (0 == count_in_grid) - mlog << Warning << method_name << " no valid point observation data within to grid\n"; + mlog << Warning << "\n" << method_name + << "no valid point observation data within to grid\n\n"; else { status = true; - mlog << Debug(3) << method_name << " count in grid: " << count_in_grid + mlog << Debug(3) << method_name << "count in grid: " << count_in_grid << " out of " << obs_count << " (" << ((obs_count > 0) ? 1.0*count_in_grid/obs_count*100 : 0) << "%)\n"; } @@ -1828,14 +1836,13 @@ static void get_grid_mapping_latlon( 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(lats, lons) -> "; int *to_cell_counts = new int[to_size]; int *mapping_indices = new int[data_size]; for (int xIdx=0; xIdx= to_size ) { - mlog << Error << method_name << "the mapped cell is out of range: " - << to_offset << " at " << xIdx << "\n"; + mlog << Error << "\n" << method_name + << "the mapped cell is out of range: " + << to_offset << " at " << xIdx << "\n\n"; + exit(1); } else cellMapping[to_offset].add(xIdx); } delete [] to_cell_counts; delete [] mapping_indices; - mlog << Debug(3) << method_name << " within grid: " << count_in_grid + 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"; @@ -1909,7 +1919,7 @@ static bool get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, DataPlane from_dp, to_dp; ConcatString cur_coord_name; clock_t start_clock = clock(); - static const char *method_name = "get_grid_mapping(var_lat, var_lon) "; + static const char *method_name = "get_grid_mapping(var_lat, var_lon) -> "; int to_lat_count = to_grid.ny(); int to_lon_count = to_grid.nx(); @@ -1918,7 +1928,7 @@ static bool get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, // 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 + mlog << Debug(4) << method_name << "data_size (ny*nx): " << data_size << " = " << from_lat_count << " * " << from_lon_count << "\n" << " target grid (nx,ny): " << to_lon_count * to_lat_count << " = " << to_lon_count << " * " << to_lat_count << "\n"; @@ -1926,10 +1936,16 @@ static bool get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, from_dp.set_size(from_lon_count, from_lat_count); to_dp.set_size(to_lon_count, to_lat_count); - if (IS_INVALID_NC(var_lat)) - mlog << Error << method_name << " Fail to get latitudes\n"; - else if (IS_INVALID_NC(var_lon)) - mlog << Error << method_name << " Fail to get longitudes\n"; + if (IS_INVALID_NC(var_lat)) { + mlog << Error << "\n" << method_name + << "Fail to get latitudes!\n\n"; + exit(1); + } + else if (IS_INVALID_NC(var_lon)) { + mlog << Error << "\n" << method_name + << "Fail to get longitudes!\n\n"; + exit(1); + } else if (data_size > 0) { float *latitudes = new float[data_size]; float *longitudes = new float[data_size]; @@ -1952,7 +1968,7 @@ static bool get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, static unixtime find_valid_time(NcVar time_var) { unixtime valid_time = bad_data_int; - static const char *method_name = "find_valid_time() "; + static const char *method_name = "find_valid_time() -> "; if( IS_VALID_NC(time_var) || get_dim_count(&time_var) < 2) { int time_count = get_dim_size(&time_var, 0); @@ -1962,14 +1978,14 @@ static unixtime find_valid_time(NcVar time_var) { valid_time = compute_unixtime(&time_var, time_values[0]); } else { - mlog << Error << "\n" << method_name << "-> " + mlog << Error << "\n" << method_name << "Can not read \"" << GET_NC_NAME(time_var) << "\" variable from \"" << InputFilename << "\"\n\n"; } } if (valid_time == bad_data_int) { - mlog << Error << "\n" << method_name << "-> " + mlog << Error << "\n" << method_name << "trouble finding time variable from \"" << InputFilename << "\"\n\n"; exit(1); @@ -2003,7 +2019,7 @@ ConcatString get_goes_grid_input(MetConfig config, Grid fr_grid, Grid to_grid) { void get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, ConcatString geostationary_file) { - static const char *method_name = "get_grid_mapping() "; + static const char *method_name = "get_grid_mapping() -> "; DataPlane from_dp, to_dp; ConcatString cur_coord_name; @@ -2032,7 +2048,7 @@ void get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, // Override the from nx & ny from NetCDF if exists NcFile *coord_nc_in = (NcFile *)0; if (has_coord_input) { - mlog << Debug(2) << method_name << " Reading coord file: " << cur_coord_name << "\n"; + mlog << Debug(2) << method_name << "Reading coord file: " << cur_coord_name << "\n"; coord_nc_in = open_ncfile(cur_coord_name.c_str()); if (IS_VALID_NC_P(coord_nc_in)) { from_lat_count = get_lat_count(coord_nc_in); @@ -2040,7 +2056,7 @@ void get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, } } int data_size = from_lat_count * from_lon_count; - mlog << Debug(4) << method_name << " data_size (ny*nx): " << data_size + mlog << Debug(4) << method_name << "data_size (ny*nx): " << data_size << " = " << from_lat_count << " * " << from_lon_count << "\n" << " target grid (nx,ny): " << to_lon_count * to_lat_count << " = " << to_lon_count << " * " << to_lat_count << "\n"; @@ -2097,9 +2113,10 @@ void get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, if ((latitudes[idx] > MISSING_LATLON) && (tmp_lats[idx] > MISSING_LATLON)) { if (!is_eq(latitudes[idx], tmp_lats[idx], loose_tol)) { lat_mis_matching_count++; - mlog << Warning << method_name << "diff lat at " << idx - << " binary-computing: " << latitudes[idx] << " - " - << tmp_lats[idx] << " = " << (latitudes[idx]-tmp_lats[idx]) << "\n"; + mlog << Warning << "\n" << method_name + << "diff lat at " << idx << " binary-computing: " + << latitudes[idx] << " - " << tmp_lats[idx] << " = " + << (latitudes[idx]-tmp_lats[idx]) << "\n\n"; } else lat_matching_count++; } @@ -2107,18 +2124,20 @@ void get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, if ((longitudes[idx] > MISSING_LATLON) && (tmp_lons[idx] > MISSING_LATLON)) { if (!is_eq(longitudes[idx], tmp_lons[idx], loose_tol)) { lon_mis_matching_count++; - mlog << Warning << method_name << "diff lon at " << idx - << " binary-computing: " << longitudes[idx] << " - " - << tmp_lons[idx] << " = " << (longitudes[idx]-tmp_lons[idx]) << "\n"; + mlog << Warning << "\n" << method_name + << "diff lon at " << idx << " binary-computing: " + << longitudes[idx] << " - " << tmp_lons[idx] << " = " + << (longitudes[idx]-tmp_lons[idx]) << "\n\n"; } else lon_matching_count++; } else lon_matching_count++; } if ((lon_mis_matching_count > 0) || (lat_mis_matching_count > 0)) { - mlog << Warning << method_name << "mis-matching: lon = " << lon_mis_matching_count + mlog << Warning << "\n" << method_name + << "mis-matching: lon = " << lon_mis_matching_count << " lat = " << lat_mis_matching_count << " matched: lon = " - << lon_matching_count << " lat = " << lat_matching_count << "\n"; + << lon_matching_count << " lat = " << lat_matching_count << "\n\n"; } } } @@ -2135,10 +2154,14 @@ void get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, } } } - if (0 == latitudes) - mlog << Error << method_name << " Fail to get latitude\n"; - else if (0 == longitudes) - mlog << Error << method_name << " Fail to get longitudes\n"; + if (0 == latitudes) { + mlog << Error << "\n" << method_name + << "Fail to get latitudes!\n\n"; + } + else if (0 == longitudes) { + mlog << Error << "\n" << method_name + << "Fail to get longitudes!\n\n"; + } else { check_lat_lon(data_size, latitudes, longitudes); get_grid_mapping_latlon(from_dp, to_dp, to_grid, cellMapping, latitudes, @@ -2185,7 +2208,8 @@ static NcVar get_goes_nc_var(NcFile *nc, const ConcatString var_name, var_data = get_nc_var(nc, var_name.split("_")[0].c_str()); } if (IS_INVALID_NC(var_data)) { - mlog << Error << "The variable \"" << var_name << "\" does not exist\n"; + mlog << Error << "\nget_goes_nc_var() -> " + << "The variable \"" << var_name << "\" does not exist\n\n"; if (exit_if_error) exit(1); } return var_data; @@ -2245,7 +2269,10 @@ static bool is_time_mismatch(NcFile *nc_in, NcFile *nc_adp) { << aod_end_time << " != " << adp_end_time << "\n"; } if (time_mismatch) { - mlog << Error << "The ADP input does not match with AOD input\n" << log_message; + mlog << Error << "\nis_time_mismatch() -> " + << "The ADP input does not match with AOD input: " + << log_message << "\n\n"; + exit(1); } return time_mismatch; } @@ -2272,7 +2299,7 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, uchar *adp_data = new uchar[from_data_size]; float *from_data = new float[from_data_size]; unsigned short *adp_qc_data = new unsigned short[from_data_size]; - static const char *method_name = "regrid_goes_variable() "; + static const char *method_name = "regrid_goes_variable() -> "; // -99 is arbitrary number as invalid QC value memset(qc_data, -99, from_data_size*sizeof(uchar)); @@ -2281,10 +2308,6 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, NcVar var_adp; NcVar var_adp_qc; NcVar var_data = get_goes_nc_var(nc_in, vinfo->name()); - //if (IS_INVALID_NC(var_data)) { - // mlog << Error << "The variable \"" << vinfo->name() << "\" does not exist\n"; - // exit(1); - //} bool is_dust_only = false; bool is_smoke_only = false; string actual_var_name = GET_NC_NAME(var_data); @@ -2315,11 +2338,15 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, << " for " << GET_NC_NAME(var_adp) << ".\n"; } } - else mlog << Warning << method_name << "found QC var name (" << qc_var_name - << " for " << GET_NC_NAME(var_adp) << ") but does not exist.\n"; + else { + mlog << Warning << "\n" << method_name + << "QC var name (" << qc_var_name + << " for " << GET_NC_NAME(var_adp) + << ") does not exist.\n\n"; + } } } - mlog << Debug(5) << method_name << " is_dust: " << is_dust_only + mlog << Debug(5) << method_name << "is_dust: " << is_dust_only << ", is_smoke: " << is_smoke_only << "\n"; //AOD:ancillary_variables = "DQF" ; byte DQF(y, x) ; @@ -2330,8 +2357,11 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, has_qc_var = true; mlog << Debug(3) << method_name << "found QC var: " << qc_var_name << ".\n"; } - else mlog << Warning << method_name << "found QC var name (" - << qc_var_name << ") but does not exist.\n"; + else { + mlog << Warning << "\n" << method_name + << "QC var name (" << qc_var_name + << ") does not exist.\n"; + } } get_nc_data(&var_data, (float *)from_data); @@ -2352,7 +2382,6 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, int non_missing_count = 0; int qc_filtered_count = 0; int adp_qc_filtered_count = 0; - // int global_attr_count; float data_value; float from_min_value = 10e10; float from_max_value = -10e10; @@ -2388,7 +2417,7 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, if (from_max_value < data_value) from_max_value = data_value; } - //Filter by QC flag + // Filter by QC flag qc_value = qc_data[from_index]; if (!has_qc_var || !has_qc_flags || qc_flags.has(qc_value)) { for(int i=0; icensor_thresh().n(); i++) { @@ -2438,7 +2467,7 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, if (0 == data_count % 2) to_value = (to_value + dataArray[(data_count/2)+1])/2; } - else to_value = dataArray.sum() / data_count; // UW_Mean + else to_value = dataArray.sum() / data_count; // UW_Mean to_dp.set(to_value, xIdx, yIdx); to_cell_count++; @@ -2459,7 +2488,7 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, delete [] from_data; delete [] adp_qc_data; - mlog << Debug(4) << method_name << " Count: actual: " << to_cell_count + mlog << Debug(4) << method_name << "Count: actual: " << to_cell_count << ", missing: " << missing_count << ", non_missing: " << non_missing_count << "\n\tFiltered: by QC: " << qc_filtered_count << ", by adp QC: " << adp_qc_filtered_count @@ -2468,8 +2497,10 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, << "\n\tRange: data: [" << from_min_value << " - " << from_max_value << "] QC: [" << qc_min_value << " - " << qc_max_value << "]\n"; - if (to_cell_count == 0) - mlog << Warning << method_name << " No valid data\n"; + if (to_cell_count == 0) { + mlog << Warning << "\n" << method_name + << "No valid data!\n\n"; + } mlog << Debug(LEVEL_FOR_PERFORMANCE) << method_name << "took " << (clock()-start_clock)/double(CLOCKS_PER_SEC) << " seconds\n"; @@ -2483,7 +2514,7 @@ static void save_geostationary_data(const ConcatString geostationary_file, bool has_error = false; int deflate_level = 0; clock_t start_clock = clock(); - static const char *method_name = "save_geostationary_data() "; + static const char *method_name = "save_geostationary_data() -> "; NcFile *nc_file = open_ncfile(geostationary_file.text(), true); NcDim xdim = add_dim(nc_file, dim_name_lon, grid_data.nx); @@ -2504,7 +2535,8 @@ static void save_geostationary_data(const ConcatString geostationary_file, add_att(&lat_var, "dy_rad", grid_data.dy_rad); if(!put_nc_data((NcVar *)&lat_var, latitudes)) { has_error = true; - mlog << Error << "Can not save latitudes\n"; + mlog << Warning << "\nsave_geostationary_data() -> " + << "Cannot save latitudes!\n\n"; } } if (IS_VALID_NC(lon_var)) { @@ -2519,7 +2551,8 @@ static void save_geostationary_data(const ConcatString geostationary_file, add_att(&lon_var, "dx_rad", grid_data.dx_rad); if(!put_nc_data((NcVar *)&lon_var, longitudes)) { has_error = true; - mlog << Error << "Can not save longitudes\n"; + mlog << Warning << "\nsave_geostationary_data() -> " + << "Cannot save longitudes!\n\n"; } } @@ -2527,7 +2560,8 @@ static void save_geostationary_data(const ConcatString geostationary_file, if (has_error) { remove(geostationary_file.c_str()); - mlog << Warning << "The geostationary data file (" + mlog << Warning << "\nsave_geostationary_data() -> " + << "The geostationary data file (" << geostationary_file << ") was not saved!\n"; } else { @@ -2578,7 +2612,7 @@ bool has_lat_lon_vars(NcFile *nc) { bool has_lon_var = IS_VALID_NC(get_nc_var_lon(nc)); bool has_time_var = IS_VALID_NC(get_nc_var_time(nc)); - //TODO: chech if this is a gridded data or a point data here!!! + //TODO: check if this is a gridded data or a point data here!!! mlog << Debug(7) << "has_lat_lon_vars() " << " has_lat_var: " << has_lat_var << ", has_lon_var: " << has_lon_var diff --git a/met/src/tools/tc_utils/tc_gen/tc_gen.cc b/met/src/tools/tc_utils/tc_gen/tc_gen.cc index 4b5acaf56e..a54d765945 100644 --- a/met/src/tools/tc_utils/tc_gen/tc_gen.cc +++ b/met/src/tools/tc_utils/tc_gen/tc_gen.cc @@ -467,6 +467,7 @@ void do_genesis_ctc(const TCGenVxOpt &vx_opt, diff.DevCategory = FNOYGenesis; diff.OpsCategory = FNOYGenesis; } + // Matched genesis pairs (DISCARD, HIT, or FALSE ALARM) else { @@ -690,14 +691,13 @@ void get_atcf_files(const StringArray &source, void process_fcst_tracks(const StringArray &files, const StringArray &model_suffix, GenesisInfoArray &fcst_ga) { - int i, j, k; + int i, j; int n_lines, tot_lines, tot_tracks, n_genesis; ConcatString suffix; LineDataFile f; ATCFTrackLine line; TrackInfoArray fcst_ta; GenesisInfo fcst_gi; - bool keep; int valid_freq_sec = conf_info.ValidFreqHr*sec_per_hour; @@ -803,8 +803,118 @@ void process_fcst_tracks(const StringArray &files, << fcst_ga.n() << "] " << fcst_ga[i].serialize() << "\n"; } + + // Set metadata pointer + suffix = model_suffix[i]; + line.set_tech_suffix(&suffix); + + // Process the input track lines + while(f >> line) { + + // Skip off-hour track points + if((line.valid_hour() % valid_freq_sec) != 0) continue; + + // Increment the line counter + n_lines++; + + // Store all BEST track lines + if(line.is_best_track()) { + best_ta.add(line, false, true); + } + // Store only 0-hour operational track lines + else if(line.is_oper_track() && line.lead() == 0) { + oper_ta.add(line); + } + } + + // Close the current file + f.close(); + + } // end for i + + // Dump out the total number of lines + mlog << Debug(3) + << "Found a total of " << best_ta.n() << " " + << conf_info.BestEventInfo.Technique + << " tracks and " << oper_ta.n() << " " + << conf_info.OperTechnique + << " operational tracks, from " << n_lines + << " track lines, from " << files.n() + << " input files.\n"; + + // Dump out very verbose output + if(mlog.verbosity_level() >= 6) { + mlog << Debug(6) << "BEST tracks:\n" + << best_ta.serialize_r() << "\n" + << "Operational tracks:\n" + << oper_ta.serialize_r() << "\n"; } + // Search the BEST tracks for genesis events + for(i=0; i " + << case_cs << "neither " << best_ga[i_bga].storm_id() + << " nor " << best_gi.storm_id() + << " matches the basin!\n\n"; + continue; + } + } + + // Compute the distance to land + best_gi.set_dland(conf_info.compute_dland( + best_gi.lat(), -1.0*best_gi.lon())); + + // Store the genesis event + best_ga.add(best_gi); + + } // end for i + + // Dump out the number of genesis events + mlog << Debug(2) << "Found " << best_ga.n() + << " BEST genesis events.\n"; + return; } @@ -1328,7 +1438,8 @@ void write_genmpr_row(StatHdrColumns &shc, void write_nc(GenCTCInfo &gci) { int i; ConcatString var_name, long_name; - unixtime valid_beg, valid_end; + unixtime valid_beg = (unixtime) 0; + unixtime valid_end = (unixtime) 0; // Allocate memory float *data = (float *) 0; diff --git a/met/src/tools/tc_utils/tc_gen/tc_gen_conf_info.cc b/met/src/tools/tc_utils/tc_gen/tc_gen_conf_info.cc index 13bcab8455..cdab8c11da 100644 --- a/met/src/tools/tc_utils/tc_gen/tc_gen_conf_info.cc +++ b/met/src/tools/tc_utils/tc_gen/tc_gen_conf_info.cc @@ -580,7 +580,6 @@ void TCGenConfInfo::read_config(const char *default_file_name, void TCGenConfInfo::process_config() { Dictionary *dict = (Dictionary *) 0; TCGenVxOpt vx_opt; - bool status; int i, beg, end; // Conf: init_freq @@ -772,7 +771,8 @@ double TCGenConfInfo::compute_dland(double lat, double lon) { //////////////////////////////////////////////////////////////////////// ConcatString TCGenConfInfo::compute_basin(double lat, double lon) { - double x_dbl, y_dbl, dist; + double x_dbl, y_dbl; + int x, y, i; // Load the basin data, if needed.