From 45a16ef0d239782ad6097f045a9323cf13516261 Mon Sep 17 00:00:00 2001 From: John Halley Gotway Date: Mon, 1 Feb 2021 12:25:23 -0700 Subject: [PATCH 1/6] Per #1644, no actual code changes here. Just formatting and spacing. For example, replace double ;; with single ;' --- met/src/tools/other/point2grid/point2grid.cc | 20 ++++++++++++-------- 1 file changed, 12 insertions(+), 8 deletions(-) diff --git a/met/src/tools/other/point2grid/point2grid.cc b/met/src/tools/other/point2grid/point2grid.cc index 70f1c7b117..f24c755e57 100644 --- a/met/src/tools/other/point2grid/point2grid.cc +++ b/met/src/tools/other/point2grid/point2grid.cc @@ -405,10 +405,12 @@ void process_data_file() { // Open the output file open_nc(to_grid, run_cs); - if (goes_data) + if (goes_data) { process_goes_file(nc_in, config, vinfo, fr_grid, to_grid); - else if (TYPE_OBS == obs_type) + } + else if (TYPE_OBS == obs_type) { process_point_file(nc_in, config, vinfo, fr_grid, to_grid); + } else if (TYPE_NCCF == obs_type) { process_point_nccf_file(nc_in, config, vinfo, fr_mtddf, to_grid); unsetenv(nc_att_met_point_nccf); @@ -723,7 +725,9 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, obs_count_zero_to = obs_count_non_zero_to = 0; obs_count_zero_from = obs_count_non_zero_from = 0; for(int i=0; iclear(); @@ -743,7 +747,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, else { if (use_var_id) { if (!var_names.has(vname, var_idx_or_gc)) { - exit_by_field_name_error = true;; + exit_by_field_name_error = true; error_msg << "The variable \"" << vname << "\" is not available.\n"; } } @@ -758,7 +762,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, sprintf(grib_code, "%d", var_idx_or_gc); } else { - exit_by_field_name_error = true;; + exit_by_field_name_error = true; error_msg << "Invalid GRIB code [" << vname << "]\n"; } } @@ -771,7 +775,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, } } if (not_found_grib_code) { - exit_by_field_name_error = true;; + exit_by_field_name_error = true; error_msg << "No data for the GRIB code [" << vname << "]\n"; } } @@ -877,7 +881,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, continue; } - //Filter by QC flag + // Filter by QC flag if (has_qc_flags && !qc_idx_array.has(obs_qty_ids[idx])) { filtered_by_qc++; continue; @@ -1923,7 +1927,7 @@ static bool get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, int to_lat_count = to_grid.ny(); int to_lon_count = to_grid.nx(); - int from_lat_count = fr_grid.ny();; + int from_lat_count = fr_grid.ny(); int from_lon_count = fr_grid.nx(); // Override the from nx & ny from NetCDF if exists @@ -2026,7 +2030,7 @@ void get_grid_mapping(Grid fr_grid, Grid to_grid, IntArray *cellMapping, clock_t start_clock = clock(); int to_lat_count = to_grid.ny(); int to_lon_count = to_grid.nx(); - int from_lat_count = fr_grid.ny();; + int from_lat_count = fr_grid.ny(); int from_lon_count = fr_grid.nx(); bool has_coord_input = false; From b260607ddb9455bbaa13d1684f2c1924524c6f8a Mon Sep 17 00:00:00 2001 From: John Halley Gotway Date: Mon, 1 Feb 2021 12:27:05 -0700 Subject: [PATCH 2/6] Per #1644, FOUND THE BUG! It's a copy/paste error. We had var_name_map.end() that should be def_var_name_map.end(). Fixing that gets rid of the runtime hang.' --- met/src/tools/other/point2grid/point2grid_conf_info.cc | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/met/src/tools/other/point2grid/point2grid_conf_info.cc b/met/src/tools/other/point2grid/point2grid_conf_info.cc index 2dd1bc1245..73c1f3a78b 100644 --- a/met/src/tools/other/point2grid/point2grid_conf_info.cc +++ b/met/src/tools/other/point2grid/point2grid_conf_info.cc @@ -50,7 +50,7 @@ void PointToGridConfInfo::clear() { quality_mark_thresh = bad_data_int; version.clear(); valid_time = 0; - //def_var_name_map.clear(); + def_var_name_map.clear(); var_name_map.clear(); return; @@ -126,7 +126,7 @@ void PointToGridConfInfo::process_config() { ConcatString PointToGridConfInfo::get_var_id(const ConcatString var_name) { ConcatString var_id; - + map::iterator ptr; for (ptr=var_name_map.begin(); ptr != var_name_map.end(); ptr++) { if( ptr->second == var_name ) { @@ -134,14 +134,16 @@ ConcatString PointToGridConfInfo::get_var_id(const ConcatString var_name) { break; } } + if( var_id.empty() ) { - for (ptr=def_var_name_map.begin(); ptr != var_name_map.end(); ptr++) { + for (ptr=def_var_name_map.begin(); ptr != def_var_name_map.end(); ptr++) { if( ptr->second == var_name ) { var_id = ptr->first; break; } } } + return var_id; } From 4c3942769b91e61b7e4c31fbb1d5cc94fea33edc Mon Sep 17 00:00:00 2001 From: Howard Soh Date: Mon, 8 Feb 2021 20:58:29 -0700 Subject: [PATCH 3/6] #1630 Display a warning instead of error message with invalid variable if the input data is empty --- met/src/tools/other/point2grid/point2grid.cc | 67 +++++++++++--------- 1 file changed, 38 insertions(+), 29 deletions(-) diff --git a/met/src/tools/other/point2grid/point2grid.cc b/met/src/tools/other/point2grid/point2grid.cc index f24c755e57..f5a5f57e3b 100644 --- a/met/src/tools/other/point2grid/point2grid.cc +++ b/met/src/tools/other/point2grid/point2grid.cc @@ -659,6 +659,7 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, // Check and read obs_vid and obs_var if exists bool success_to_read = true; + bool empty_input = (nhdr == 0 && nobs == 0); NcVar obs_vid_var = get_var(nc_in, nc_var_obs_vid); if (IS_VALID_NC(obs_vid_var)) { if (success_to_read) success_to_read = get_nc_data_int_array(nc_in, nc_var_obs_vid, obs_ids); @@ -673,30 +674,32 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, } } - if (success_to_read) - success_to_read = get_nc_data_int_array(nc_in, nc_var_obs_hid, obs_hids); - if (success_to_read) - success_to_read = get_nc_data_int_array(nc_in, nc_var_hdr_vld, hdr_vld_ids); - if (success_to_read) - success_to_read = get_nc_data_int_array(nc_in, nc_var_hdr_typ, hdr_typ_ids); - if (success_to_read) - success_to_read = get_nc_data_int_array(nc_in, nc_var_obs_qty, obs_qty_ids); - if (success_to_read) - 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 (success_to_read) - success_to_read = get_nc_data_float_array(nc_in, nc_var_obs_lvl, obs_lvls); - if (success_to_read) - success_to_read = get_nc_data_float_array(nc_in, nc_var_obs_hgt, obs_hgts); - if (success_to_read) - success_to_read = get_nc_data_float_array(nc_in, nc_var_obs_val, obs_vals); - if (success_to_read) - success_to_read = get_nc_data_string_array( - nc_in, nc_var_hdr_vld_tbl, &hdr_valid_times); - if (success_to_read) - success_to_read = get_nc_data_string_array( - nc_in, nc_var_obs_qty_tbl, &qc_tables); + if (!empty_input) { + if (success_to_read) + success_to_read = get_nc_data_int_array(nc_in, nc_var_obs_hid, obs_hids); + if (success_to_read) + success_to_read = get_nc_data_int_array(nc_in, nc_var_hdr_vld, hdr_vld_ids); + if (success_to_read) + success_to_read = get_nc_data_int_array(nc_in, nc_var_hdr_typ, hdr_typ_ids); + if (success_to_read) + success_to_read = get_nc_data_int_array(nc_in, nc_var_obs_qty, obs_qty_ids); + if (success_to_read) + 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 (success_to_read) + success_to_read = get_nc_data_float_array(nc_in, nc_var_obs_lvl, obs_lvls); + if (success_to_read) + success_to_read = get_nc_data_float_array(nc_in, nc_var_obs_hgt, obs_hgts); + if (success_to_read) + success_to_read = get_nc_data_float_array(nc_in, nc_var_obs_val, obs_vals); + if (success_to_read) + success_to_read = get_nc_data_string_array( + nc_in, nc_var_hdr_vld_tbl, &hdr_valid_times); + if (success_to_read) + success_to_read = get_nc_data_string_array( + nc_in, nc_var_obs_qty_tbl, &qc_tables); + } if (success_to_read) { bool has_qc_flags = (qc_flags.n() > 0); IntArray qc_idx_array = prepare_qc_array(qc_flags, qc_tables); @@ -801,11 +804,17 @@ 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"; - exit(1); + if (empty_input) { + mlog << Warning << "\n" << method_name + << error_msg << "\tBut ignored because of empty input\n\n"; + } + else { + 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"; + exit(1); + } } // Check the time range. Apply the time window From 530a874adc2668de736fd134fbb9232858b55d9a Mon Sep 17 00:00:00 2001 From: johnhg Date: Wed, 10 Feb 2021 10:51:49 -0700 Subject: [PATCH 4/6] Feature 1658 grib_tables (#1659) * Per #1658, update MXUPHL entries. * Per #1658, updating long name for MAXREF, MAXUVV, and MAXDVV. --- met/data/table_files/grib1_ncep_129_7.txt | 8 ++++---- met/data/table_files/grib2_all.txt | 10 +++++----- 2 files changed, 9 insertions(+), 9 deletions(-) diff --git a/met/data/table_files/grib1_ncep_129_7.txt b/met/data/table_files/grib1_ncep_129_7.txt index aec86c19f3..eef20ff595 100644 --- a/met/data/table_files/grib1_ncep_129_7.txt +++ b/met/data/table_files/grib1_ncep_129_7.txt @@ -234,10 +234,10 @@ GRIB1 232 129 7 -1 "OZMAX8" "Ozone Daily Max from 8-hour Average" "ppbV" 233 129 7 -1 "PDMAX1" "PM 2.5 Daily Max from 1-hour Average" "mcg/m^3" 234 129 7 -1 "PDMX24" "PM 2.5 Daily Max from 24-hour Average" "mcg/m^3" -235 129 7 -1 "MAXREF" "Hourly Maximum of Simulated Reflectivity at 1 km AGL" "dbZ" -236 129 7 -1 "MXUPHL" "Hourly Maximum of Updraft Helicity over layer 2km to 5 km AGL" "m^2/s^2" -237 129 7 -1 "MAXUVV" "Hourly Maximum of Upward Vertical Velocity in the lowest 400hPa" "m/s" -238 129 7 -1 "MAXDVV" "Hourly Maximum of Downward Vertical Velocity in the lowest 400hPa" "m/s" +235 129 7 -1 "MAXREF" "Hourly Maximum of Simulated Reflectivity" "dbZ" +236 129 7 -1 "MXUPHL" "Hourly Maximum of Updraft Helicity" "m^2/s^2" +237 129 7 -1 "MAXUVV" "Hourly Maximum of Upward Vertical Velocity" "m/s" +238 129 7 -1 "MAXDVV" "Hourly Maximum of Downward Vertical Velocity" "m/s" 239 129 7 -1 "MAXVIG" "Hourly Maximum of Column Vertical Integrated Graupel" "kg/m^2" 240 129 7 -1 "RETOP" "Radar Echo Top (18.3 DBZ)" "m" 241 129 7 -1 "VRATE" "Ventilation Rate" "m^2/s" diff --git a/met/data/table_files/grib2_all.txt b/met/data/table_files/grib2_all.txt index c41a244178..270675890f 100644 --- a/met/data/table_files/grib2_all.txt +++ b/met/data/table_files/grib2_all.txt @@ -421,11 +421,11 @@ GRIB2 0 1 0 255 0 0 1 97 "MASSDS" "Mass Density of Snow" "kg/m^3" 0 1 0 255 0 0 20 2 "MASSMR" "Mass Mixing Ratio (Mass Fraction in Air)" "kg/kg" 0 1 0 255 0 0 1 28 "MAXAH" "Maximum Absolute Humidity" "kg/m^3" -0 0 0 255 7 1 2 221 "MAXDVV" "Hourly Maximum of Downward Vertical Velocity in the lowest 400hPa" "m/s" +0 0 0 255 7 1 2 221 "MAXDVV" "Hourly Maximum of Downward Vertical Velocity" "m/s" 0 1 0 255 0 0 2 21 "MAXGUST" "Maximum Wind Speed" "m/s" -0 0 0 255 7 1 16 198 "MAXREF" "Hourly Maximum of Simulated Reflectivity at 1 km AGL" "dB" +0 0 0 255 7 1 16 198 "MAXREF" "Hourly Maximum of Simulated Reflectivity" "dB" 0 1 0 255 0 0 1 27 "MAXRH" "Maximum Relative Humidity" "%" -0 0 0 255 7 1 2 220 "MAXUVV" "Hourly Maximum of Upward Vertical Velocity in the lowest 400hPa" "m/s" +0 0 0 255 7 1 2 220 "MAXUVV" "Hourly Maximum of Upward Vertical Velocity" "m/s" 0 0 0 255 7 1 2 222 "MAXUW" "U Component of Hourly Maximum 10m Wind Speed" "m/s" 0 0 0 255 7 1 2 223 "MAXVW" "V Component of Hourly Maximum 10m Wind Speed" "m/s" 10 1 0 255 0 0 0 24 "MAXWH" "Maximum Individual Wave Height" "m" @@ -460,7 +460,7 @@ GRIB2 0 1 0 255 0 0 19 28 "MWTURB" "Mountain Wave Turbulence (Eddy Dissipation Rate)" "m2/3s-1" 0 0 0 255 7 1 19 192 "MXSALB" "Maximum Snow Albedo" "%" 0 1 0 255 0 0 19 17 "MXSALB" "Maximum Snow Albedosee Note 1" "%" -0 0 0 255 7 1 7 199 "MXUPHL" "Hourly Maximum of Updraft Helicity over Layer 2km to 5 km AGL" "m^2/s^2" +0 0 0 255 7 1 7 199 "MXUPHL" "Hourly Maximum of Updraft Helicity" "m^2/s^2" 10 1 0 255 0 0 0 30 "MZPTSW" "Mean Zero-Crossing Period of The Total Swell" "s" 10 1 0 255 0 0 0 29 "MZPWW" "Mean Zero-Crossing Period of The Wind Waves" "s" 10 1 0 255 0 0 0 28 "MZWPER" "Mean Zero-Crossing Wave Period" "s" @@ -1016,4 +1016,4 @@ GRIB2 4 1 0 255 0 0 6 1 "XLONG" "Solar X-ray Flux (XRS Long)" "W/m^2" 4 1 0 255 0 0 8 0 "XRAYRAD" "X-Ray Radiance" "W/sr/m^2" 4 1 0 255 0 0 6 2 "XSHRT" "Solar X-ray Flux (XRS Short)" "W/m^2" -10 1 0 255 0 0 2 10 "ZVCICEP" "Zonal Vector Component of Vertically Integrated Ice Internal Pressure" "Pa*m" \ No newline at end of file +10 1 0 255 0 0 2 10 "ZVCICEP" "Zonal Vector Component of Vertically Integrated Ice Internal Pressure" "Pa*m" From 4c9ecdee283ea0ed2b1c0bbccf2d443962d12c51 Mon Sep 17 00:00:00 2001 From: "Julie.Prestopnik" Date: Thu, 11 Feb 2021 15:22:50 -0700 Subject: [PATCH 5/6] Modified format of release notes --- met/docs/Users_Guide/release-notes.rst | 22 +++++++++++++--------- 1 file changed, 13 insertions(+), 9 deletions(-) diff --git a/met/docs/Users_Guide/release-notes.rst b/met/docs/Users_Guide/release-notes.rst index 4c3eb561bb..679c9e7fbf 100644 --- a/met/docs/Users_Guide/release-notes.rst +++ b/met/docs/Users_Guide/release-notes.rst @@ -29,22 +29,26 @@ Version `10.0.0-beta3 `_ release no * When reading MET NetCDF files, parse the "init_time" and "valid_time" attributes (`#1346 `_). * Python embedding enhancements: - * Complete support for Python XArray embedding (`#1534 `_). - * Correct error messages from Python embedding (`#1473 `_). + + * Complete support for Python XArray embedding (`#1534 `_). + * Correct error messages from Python embedding (`#1473 `_). * Application code: * Enhance Plot-Point-Obs to support regridding in the config file (`#1627 `_). * Update ASCII2NC and Point2Grid to create empty output files for zero input observations instead of erroring out (`#1630 `_). * Point2Grid Tool: - * Improve the Point2Grid runtime performance (`#1421 `_). - * Process point observations by variable name instead of GRIB code (`#1408 `_). - * Support the 2-dimensional time variable in Himawari data files (`#1580 `_). + + * Improve the Point2Grid runtime performance (`#1421 `_). + * Process point observations by variable name instead of GRIB code (`#1408 `_). + * Support the 2-dimensional time variable in Himawari data files (`#1580 `_). + * TC-Gen Tool: - * Overhaul the Tropical Cyclone genesis matching logic, add the development and operational scoring algorithms, and add many config file options (`#1448 `_). - * Add config file options to filter data by initialization time (init_inc and init_exc) and hurricane basin (basin_mask) (`#1626 `_). - * Add the genesis matched pair (GENMPR) output line type (`#1597 `_). - * Add a gridded NetCDF output file with counts for genesis events and track points (`#1430 `_). + + * Overhaul the Tropical Cyclone genesis matching logic, add the development and operational scoring algorithms, and add many config file options (`#1448 `_). + * Add config file options to filter data by initialization time (init_inc and init_exc) and hurricane basin (basin_mask) (`#1626 `_). + * Add the genesis matched pair (GENMPR) output line type (`#1597 `_). + * Add a gridded NetCDF output file with counts for genesis events and track points (`#1430 `_). Version `10.0.0-beta2 `_ release notes (20201207) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ From 11a5419931e3488e26e5e9359fc43efc46444478 Mon Sep 17 00:00:00 2001 From: johnhg Date: Wed, 17 Feb 2021 12:08:59 -0700 Subject: [PATCH 6/6] Feature 1450 hersbach (#1662) * Per #1450, add new ECNT columns for Hersback CRPS. Still need to actually compute the stats though. * Per #1450, update NumArray functions to only sort if the data is not yet sorted. And check for bad data when computing the standard deviation. * Per #1450, add code to compute the empirical CRPS value. * Per #1450, large change to the new output for the empirical CRPS. In order to aggregate decomposed empirical CRPS reliability and potential correctly, we'd need to write (n+1)*2 additional columns. While the empirical crps can be aggregated as a weighted mean, the decomposition cannot. It just isn't feasible to do this in the ECNT line type. If this reliability and potential really are required, recommend that we add an entirely new CRPS line type instead of tacking onto ECNT. These changes simply remove reliabilit and potential from the output. * Per #1450 and #1451, replacing single CRPS_CLIMO column with CRPSCL and CRPSCL_EMP which will be needed for #1451. * Per #1450, delete temp files I'd accidentally committed. * Per #1450, update the user's guide with CRPS updates. * Fix bug replacing crpss_emp with crpss_gaus. --- met/.cproject | 362 ------------------ met/.project | 78 ---- .../table_files/met_header_columns_V10.0.txt | 2 +- met/docs/Users_Guide/appendixC.rst | 8 +- met/docs/Users_Guide/ensemble-stat.rst | 18 +- met/docs/Users_Guide/point-stat.rst | 2 +- met/docs/Users_Guide/refs.rst | 6 + met/src/basic/vx_util/num_array.cc | 4 +- met/src/basic/vx_util/stat_column_defs.h | 4 +- met/src/libcode/vx_stat_out/stat_columns.cc | 27 +- .../libcode/vx_statistics/compute_stats.cc | 19 +- met/src/libcode/vx_statistics/ens_stats.cc | 50 ++- met/src/libcode/vx_statistics/ens_stats.h | 5 +- .../vx_statistics/pair_data_ensemble.cc | 116 ++++-- .../vx_statistics/pair_data_ensemble.h | 19 +- .../core/stat_analysis/aggr_stat_line.cc | 54 ++- .../tools/core/stat_analysis/aggr_stat_line.h | 2 +- .../core/stat_analysis/parse_stat_line.cc | 11 +- .../core/stat_analysis/parse_stat_line.h | 5 +- test/hdr/met_10_0.hdr | 2 +- 20 files changed, 247 insertions(+), 547 deletions(-) delete mode 100644 met/.cproject delete mode 100644 met/.project diff --git a/met/.cproject b/met/.cproject deleted file mode 100644 index a4fe10144f..0000000000 --- a/met/.cproject +++ /dev/null @@ -1,362 +0,0 @@ - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - diff --git a/met/.project b/met/.project deleted file mode 100644 index fcf9456036..0000000000 --- a/met/.project +++ /dev/null @@ -1,78 +0,0 @@ - - - MET_svn - - - - - - org.eclipse.cdt.managedbuilder.core.genmakebuilder - clean,full,incremental, - - - ?name? - - - - org.eclipse.cdt.make.core.append_environment - true - - - org.eclipse.cdt.make.core.autoBuildTarget - all - - - org.eclipse.cdt.make.core.buildArguments - - - - org.eclipse.cdt.make.core.buildCommand - make - - - org.eclipse.cdt.make.core.cleanBuildTarget - clean - - - org.eclipse.cdt.make.core.contents - org.eclipse.cdt.make.core.activeConfigSettings - - - org.eclipse.cdt.make.core.enableAutoBuild - false - - - org.eclipse.cdt.make.core.enableCleanBuild - true - - - org.eclipse.cdt.make.core.enableFullBuild - true - - - org.eclipse.cdt.make.core.fullBuildTarget - all - - - org.eclipse.cdt.make.core.stopOnError - true - - - org.eclipse.cdt.make.core.useDefaultBuildCmd - true - - - - - org.eclipse.cdt.managedbuilder.core.ScannerConfigBuilder - - - - - - org.eclipse.cdt.core.cnature - org.eclipse.cdt.core.ccnature - org.eclipse.cdt.managedbuilder.core.managedBuildNature - org.eclipse.cdt.managedbuilder.core.ScannerConfigNature - - diff --git a/met/data/table_files/met_header_columns_V10.0.txt b/met/data/table_files/met_header_columns_V10.0.txt index 851ae5e3ee..d2c8ae5c4f 100644 --- a/met/data/table_files/met_header_columns_V10.0.txt +++ b/met/data/table_files/met_header_columns_V10.0.txt @@ -17,7 +17,7 @@ V10.0 : STAT : PJC : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID V10.0 : STAT : PRC : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL (N_THRESH) THRESH_[0-9]* PODY_[0-9]* POFD_[0-9]* V10.0 : STAT : PSTD : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL (N_THRESH) BASER BASER_NCL BASER_NCU RELIABILITY RESOLUTION UNCERTAINTY ROC_AUC BRIER BRIER_NCL BRIER_NCU BRIERCL BRIERCL_NCL BRIERCL_NCU BSS BSS_SMPL THRESH_[0-9]* V10.0 : STAT : ECLV : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL BASER VALUE_BASER (N_PTS) CL_[0-9]* VALUE_[0-9]* -V10.0 : STAT : ECNT : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL N_ENS CRPS CRPSS IGN ME RMSE SPREAD ME_OERR RMSE_OERR SPREAD_OERR SPREAD_PLUS_OERR +V10.0 : STAT : ECNT : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL N_ENS CRPS CRPSS IGN ME RMSE SPREAD ME_OERR RMSE_OERR SPREAD_OERR SPREAD_PLUS_OERR CRPSCL CRPS_EMP CRPSCL_EMP CRPSS_EMP V10.0 : STAT : RPS : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL N_PROB RPS_REL RPS_RES RPS_UNC RPS RPSS RPSS_SMPL RPS_COMP V10.0 : STAT : RHIST : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL (N_RANK) RANK_[0-9]* V10.0 : STAT : PHIST : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL BIN_SIZE (N_BIN) BIN_[0-9]* diff --git a/met/docs/Users_Guide/appendixC.rst b/met/docs/Users_Guide/appendixC.rst index d5ad94b01a..05905e9c67 100644 --- a/met/docs/Users_Guide/appendixC.rst +++ b/met/docs/Users_Guide/appendixC.rst @@ -887,9 +887,9 @@ ________________________________________________ CRPS ~~~~ -Called "CRPS" in ECNT output :numref:`table_ES_header_info_es_out_ECNT` +Called "CRPS", "CRPSCL", "CRPS_EMP", and "CRPSCL_EMP" in ECNT output :numref:`table_ES_header_info_es_out_ECNT` -The continuous ranked probability score (CRPS) is the integral, over all possible thresholds, of the Brier scores (:ref:`Gneiting et al., 2004 `). In MET, the CRPS calculation uses a normal distribution fit to the ensemble forecasts. In many cases, use of other distributions would be better. +The continuous ranked probability score (CRPS) is the integral, over all possible thresholds, of the Brier scores (:ref:`Gneiting et al., 2004 `). In MET, the CRPS is calculated two ways: using a normal distribution fit to the ensemble forecasts (CRPS and CRPSCL), and using the empirical ensemble distribution (CRPS_EMP and CRPSCL_EMP). In some cases, use of other distributions would be better. WARNING: The normal distribution is probably a good fit for temperature and pressure, and possibly a not horrible fit for winds. However, the normal approximation will not work on most precipitation forecasts and may fail for many other atmospheric variables. @@ -906,12 +906,14 @@ The score can be interpreted as a continuous version of the mean absolute error CRPS Skill Score ~~~~~~~~~~~~~~~~ -Called "CRPSS" in ECNT output :numref:`table_ES_header_info_es_out_ECNT` +Called "CRPSS" and "CRPSS_EMP" in ECNT output :numref:`table_ES_header_info_es_out_ECNT` The continuous ranked probability skill score (CRPSS) is similar to the MSESS and the BSS, in that it compares its namesake score to that of a reference forecast to produce a positively oriented score between 0 and 1. .. math:: \text{CRPSS} = 1 - \frac{\text{CRPS}_{fcst}}{ \text{CRPS}_{ref}} +For the normal distribution fit (CRPSS), the reference CRPS is computed using the climatological mean and standard deviation. For the empirical distribution (CRPSS_EMP), the reference CRPS is computed by sampling from the assumed normal climatological distribution defined by the mean and standard deviation. + IGN ~~~ diff --git a/met/docs/Users_Guide/ensemble-stat.rst b/met/docs/Users_Guide/ensemble-stat.rst index 60ef284b97..99532dcaf6 100644 --- a/met/docs/Users_Guide/ensemble-stat.rst +++ b/met/docs/Users_Guide/ensemble-stat.rst @@ -31,7 +31,7 @@ Often, the goal of ensemble forecasting is to reproduce the distribution of obse The relative position (RELP) is a count of the number of times each ensemble member is closest to the observation. For stochastic or randomly derived ensembles, this statistic is meaningless. For specified ensemble members, however, it can assist users in determining if any ensemble member is performing consistently better or worse than the others. -The ranked probability score (RPS) is included in the Ranked Probability Score (RPS) line type. It is the mean of the Brier scores computed from ensemble probabilities derived for each probability category threshold (prob_cat_thresh) specified in the configuration file. The continuous ranked probability score (CRPS) is the average the distance between the forecast (ensemble) cumulative distribution function and the observation cumulative distribution function. It is an analog of the Brier score, but for continuous forecast and observation fields. (:ref:`Gneiting et al., 2004 `). The CRPS statistic is included in the Ensemble Continuous Statistics (ECNT) line type, along with other statistics quantifying the ensemble spread and ensemble mean skill. +The ranked probability score (RPS) is included in the Ranked Probability Score (RPS) line type. It is the mean of the Brier scores computed from ensemble probabilities derived for each probability category threshold (prob_cat_thresh) specified in the configuration file. The continuous ranked probability score (CRPS) is the average the distance between the forecast (ensemble) cumulative distribution function and the observation cumulative distribution function. It is an analog of the Brier score, but for continuous forecast and observation fields. The CRPS statistic is computed using two methods: assuming a normal distribution defined by the ensemble mean and spread (:ref:`Gneiting et al., 2004 `) and using the empirical ensemble distribution (:ref:`Hersbach, 2000 `). The CRPS statistic is included in the Ensemble Continuous Statistics (ECNT) line type, along with other statistics quantifying the ensemble spread and ensemble mean skill. Ensemble observation error ~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -587,10 +587,10 @@ The format of the STAT and ASCII output of the Ensemble-Stat tool are described - Number of ensemble values * - 27 - CRPS - - The Continuous Ranked Probability Score + - The Continuous Ranked Probability Score (normal distribution) * - 28 - CRPSS - - The Continuous Ranked Probability Skill Score + - The Continuous Ranked Probability Skill Score (normal distribution) * - 29 - IGN - The Ignorance Score @@ -615,6 +615,18 @@ The format of the STAT and ASCII output of the Ensemble-Stat tool are described * - 36 - SPREAD_PLUS_OERR - The square root of the sum of unperturbed ensemble variance and the observation error variance + * - 37 + - CRPSCL + - Climatological Continuous Ranked Probability Score (normal distribution) + * - 38 + - CRPS_EMP + - The Continuous Ranked Probability Score (empirical distribution) + * - 39 + - CRPSCL_EMP + - Climatological Continuous Ranked Probability Score (empirical distribution) + * - 40 + - CRPSS_EMP + - The Continuous Ranked Probability Skill Score (empirical distribution) .. _table_ES_header_info_es_out_RPS: diff --git a/met/docs/Users_Guide/point-stat.rst b/met/docs/Users_Guide/point-stat.rst index 5756e40cd5..1ec6fdacb9 100644 --- a/met/docs/Users_Guide/point-stat.rst +++ b/met/docs/Users_Guide/point-stat.rst @@ -170,7 +170,7 @@ When the "prob" entry is set as a dictionary to define the field of interest, se Measures for comparison against climatology ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ -For each of the types of statistics mentioned above (categorical, continuous, and probabilistic), it is possible to calculate measures of skill relative to climatology. MET will accept a climatology file provided by the user, and will evaluate it as a reference forecast. Further, anomalies, i.e. departures from average conditions, can be calculated. As with all other statistics, the available measures will depend on the nature of the forecast. Common statistics that use a climatological reference include: the mean squared error skill score (MSESS), the Anomaly Correlation (ANOM_CORR and ANOM_CORR_UNCNTR), scalar and vector anomalies (SAL1L2 and VAL1L2), continuous ranked probability skill score (CRPSS), Brier Skill Score (BSS) (:ref:`Wilks, 2011 `; :ref:`Mason, 2004 `). +For each of the types of statistics mentioned above (categorical, continuous, and probabilistic), it is possible to calculate measures of skill relative to climatology. MET will accept a climatology file provided by the user, and will evaluate it as a reference forecast. Further, anomalies, i.e. departures from average conditions, can be calculated. As with all other statistics, the available measures will depend on the nature of the forecast. Common statistics that use a climatological reference include: the mean squared error skill score (MSESS), the Anomaly Correlation (ANOM_CORR and ANOM_CORR_UNCNTR), scalar and vector anomalies (SAL1L2 and VAL1L2), continuous ranked probability skill score (CRPSS and CRPSS_EMP), Brier Skill Score (BSS) (:ref:`Wilks, 2011 `; :ref:`Mason, 2004 `). Often, the sample climatology is used as a reference by a skill score. The sample climatology is the average over all included observations and may be transparent to the user. This is the case in most categorical skill scores. The sample climatology will probably prove more difficult to improve upon than a long term climatology, since it will be from the same locations and time periods as the forecasts. This may mask legitimate forecast skill. However, a more general climatology, perhaps covering many years, is often easier to improve upon and is less likely to mask real forecast skill. diff --git a/met/docs/Users_Guide/refs.rst b/met/docs/Users_Guide/refs.rst index 269242adb3..7e02dcb6d1 100644 --- a/met/docs/Users_Guide/refs.rst +++ b/met/docs/Users_Guide/refs.rst @@ -132,6 +132,12 @@ References | forecasts. *Monthly Weather Review*, 129, 550-560. | +.. _Hersbach-2000: + +| Hersbach, H., 2000: Decomposition of the Continuous Ranked Probability Score +| for Ensemble Prediction Systems. *Weather and Forecasting*, 15, 559-570. +| + .. _Jolliffe-2012: | Jolliffe, I.T., and D.B. Stephenson, 2012: *Forecast verification. A* diff --git a/met/src/basic/vx_util/num_array.cc b/met/src/basic/vx_util/num_array.cc index 0481f4516c..bf80325736 100644 --- a/met/src/basic/vx_util/num_array.cc +++ b/met/src/basic/vx_util/num_array.cc @@ -540,7 +540,7 @@ void NumArray::sort_array() { -sort(e, Nelements); +if(!Sorted) sort(e, Nelements); Sorted = true; @@ -784,7 +784,7 @@ double var; compute_mean_variance(mn, var); -stdev = square_root(var); +stdev = (is_bad_data(var) ? bad_data_double : square_root(var)); return; diff --git a/met/src/basic/vx_util/stat_column_defs.h b/met/src/basic/vx_util/stat_column_defs.h index 33d52d7156..d5a8cf0c1a 100644 --- a/met/src/basic/vx_util/stat_column_defs.h +++ b/met/src/basic/vx_util/stat_column_defs.h @@ -260,7 +260,9 @@ static const char * ecnt_columns [] = { "TOTAL", "N_ENS", "CRPS", "CRPSS", "IGN", "ME", "RMSE", "SPREAD", "ME_OERR", - "RMSE_OERR", "SPREAD_OERR", "SPREAD_PLUS_OERR" + "RMSE_OERR", "SPREAD_OERR", "SPREAD_PLUS_OERR", + "CRPSCL", "CRPS_EMP", "CRPSCL_EMP", + "CRPSS_EMP" }; static const char * rps_columns [] = { diff --git a/met/src/libcode/vx_stat_out/stat_columns.cc b/met/src/libcode/vx_stat_out/stat_columns.cc index 09905b4f08..570977fb17 100644 --- a/met/src/libcode/vx_stat_out/stat_columns.cc +++ b/met/src/libcode/vx_stat_out/stat_columns.cc @@ -3659,11 +3659,12 @@ void write_ecnt_cols(const ECNTInfo &ecnt_info, // // Ensemble Continuous Statistics // Dump out the ECNT line: - // TOTAL, N_ENS, - // CRPS, CRPSS, IGN, - // ME, RMSE, SPREAD, - // ME_OERR, RMSE_OERR, SPREAD_OERR, - // SPREAD_PLUS_OERR + // TOTAL, N_ENS, CRPS, + // CRPSS, IGN, ME, + // RMSE, SPREAD, ME_OERR, + // RMSE_OERR, SPREAD_OERR, SPREAD_PLUS_OERR, + // CRPSCL, CRPS_EMP, CRPSCL_EMP, + // CRPSS_EMP // at.set_entry(r, c+0, // Total Number of Pairs ecnt_info.n_pair); @@ -3672,10 +3673,10 @@ void write_ecnt_cols(const ECNTInfo &ecnt_info, ecnt_info.n_ens); at.set_entry(r, c+2, // Continuous Ranked Probability Score - ecnt_info.crps); + ecnt_info.crps_gaus); at.set_entry(r, c+3, // Continuous Ranked Probability Skill Score - ecnt_info.crpss); + ecnt_info.crpss_gaus); at.set_entry(r, c+4, // Ignorance Score ecnt_info.ign); @@ -3701,6 +3702,18 @@ void write_ecnt_cols(const ECNTInfo &ecnt_info, at.set_entry(r, c+11, // Mean of unperturbed spread plus observation error ecnt_info.spread_plus_oerr); + at.set_entry(r, c+12, // Gaussian climatological CRPS + ecnt_info.crpscl_gaus); + + at.set_entry(r, c+13, // Empirical ensemble CRPS + ecnt_info.crps_emp); + + at.set_entry(r, c+14, // Empirical climatological CRPS + ecnt_info.crpscl_emp); + + at.set_entry(r, c+15, // Empirical CRPSS + ecnt_info.crpss_emp); + return; } diff --git a/met/src/libcode/vx_statistics/compute_stats.cc b/met/src/libcode/vx_statistics/compute_stats.cc index 674c3c59d3..7e6e74f66f 100644 --- a/met/src/libcode/vx_statistics/compute_stats.cc +++ b/met/src/libcode/vx_statistics/compute_stats.cc @@ -1306,12 +1306,23 @@ void compute_ecnt_mean(const ECNTInfo *ecnt_info, int n, ecnt_mean.n_pair = na.sum(); // Compute unweighted mean for each statistic + for(i=0,na.erase(); i= all ensemble members + if(is_eq(obs_cdf, 0.0)) integral += (obs - fcst); + + return(integral); +} + +//////////////////////////////////////////////////////////////////////// +// +// Compute the Gaussian continuous ranked probability score (CRPS), +// the ignorance score (IGN), and probability integral transform (PIT) +// +//////////////////////////////////////////////////////////////////////// + +void compute_crps_gaus_ign_pit(double obs, const NumArray &ens_na, + double &crps, double &ign, double &pit) { double m, s, z; // Mean and standard deviation of the ensemble values @@ -1755,7 +1820,7 @@ void compute_crps_ign_pit(double obs, const NumArray &ens_na, z = (obs - m)/s; - // Compute CRPS + // Compute Gaussian CRPS crps = s*(z*(2.0*znorm(z) - 1) + 2.0*dnorm(z) - 1.0/sqrt(pi)); // Compute IGN @@ -1763,6 +1828,7 @@ void compute_crps_ign_pit(double obs, const NumArray &ens_na, // Compute PIT pit = normal_cdf(obs, m, s); + } return; @@ -1805,9 +1871,10 @@ PairDataEnsemble subset_climo_cdf_bin(const PairDataEnsemble &pd, // required for ensemble output line types. // // Include in subset: - // wgt_na, o_na, cmn_na, csd_na, cdf_na, v_na, r_na, crps_na, - // ign_na, pit_na, var_na, var_oerr_na, - // var_plus_oerr_na, mn_na, mn_oerr_na, e_na + // wgt_na, o_na, cmn_na, csd_na, cdf_na, v_na, r_na, + // crps_emp_na, crps_gaus_na, ign_na, pit_na, + // var_na, var_oerr_na, var_plus_oerr_na, + // mn_na, mn_oerr_na, e_na // // Exclude from subset: // sid_sa, lat_na, lon_na, x_na, y_na, vld_ta, lvl_ta, elv_ta, @@ -1820,7 +1887,8 @@ PairDataEnsemble subset_climo_cdf_bin(const PairDataEnsemble &pd, out_pd.cdf_na.add(pd.cdf_na[i]); out_pd.v_na.add(pd.v_na[i]); out_pd.r_na.add(pd.r_na[i]); - out_pd.crps_na.add(pd.crps_na[i]); + out_pd.crps_emp_na.add(pd.crps_emp_na[i]); + out_pd.crps_gaus_na.add(pd.crps_gaus_na[i]); out_pd.ign_na.add(pd.ign_na[i]); out_pd.pit_na.add(pd.pit_na[i]); out_pd.skip_ba.add(false); diff --git a/met/src/libcode/vx_statistics/pair_data_ensemble.h b/met/src/libcode/vx_statistics/pair_data_ensemble.h index 7fa90737bf..49a848dd95 100644 --- a/met/src/libcode/vx_statistics/pair_data_ensemble.h +++ b/met/src/libcode/vx_statistics/pair_data_ensemble.h @@ -78,11 +78,14 @@ class PairDataEnsemble : public PairBase { NumArray *e_na; // Ensemble values [n_ens][n_obs] NumArray v_na; // Number of valid ensemble values [n_obs] NumArray r_na; // Observation ranks [n_obs] - NumArray crps_na; // Continuous Ranked Probability Score [n_obs] + + NumArray crps_emp_na; // Empirical Continuous Ranked Probability Score [n_obs] + NumArray crps_gaus_na; // Gaussian CRPS [n_obs] + NumArray ign_na; // Ignorance Score [n_obs] NumArray pit_na; // Probability Integral Transform [n_obs] - int n_ens; // Number of ensemble members + int n_ens; // Number of ensemble members int n_pair; // Number of valid pairs, n_obs - sum(skip_ba) bool skip_const; // Skip cases where the observation and // all ensemble members are constant @@ -107,7 +110,11 @@ class PairDataEnsemble : public PairBase { double ssvar_bin_size; // Variance bin size for spread/skill SSVARInfo *ssvar_bins; // Ensemble spread/skill bin information [n_ssvar_bin] - double crpss; // Continuous ranked probability skill score + double crpscl_emp; // Empirical climatological CRPS score + double crpss_emp; // Empirical CRPS skill score + double crpscl_gaus; // Guassian climatological CRPS score + double crpss_gaus; // Guassian CRPS skill score + double me; // ME for ensemble mean double rmse; // RMSE for ensemble mean double me_oerr; // ME for mean of perturbed members @@ -290,8 +297,10 @@ class VxPairDataEnsemble { // //////////////////////////////////////////////////////////////////////// -extern void compute_crps_ign_pit(double, const NumArray &, double &, - double &, double &); +extern double compute_crps_emp(double, const NumArray &); + +extern void compute_crps_gaus_ign_pit(double, const NumArray &, + double &, double &, double &); // Subset pairs for a specific climatology CDF bin extern PairDataEnsemble subset_climo_cdf_bin(const PairDataEnsemble &, diff --git a/met/src/tools/core/stat_analysis/aggr_stat_line.cc b/met/src/tools/core/stat_analysis/aggr_stat_line.cc index 686a7cd87a..034afdb0cb 100644 --- a/met/src/tools/core/stat_analysis/aggr_stat_line.cc +++ b/met/src/tools/core/stat_analysis/aggr_stat_line.cc @@ -2514,7 +2514,7 @@ void aggr_ecnt_lines(LineDataFile &f, STATAnalysisJob &job, AggrENSInfo aggr; ECNTData cur; ConcatString key; - double crps_fcst, crps_climo, v; + double crps_emp, crpscl_emp, crps_gaus, crpscl_gaus, v; map::iterator it; // @@ -2553,7 +2553,8 @@ void aggr_ecnt_lines(LineDataFile &f, STATAnalysisJob &job, // if(m.count(key) == 0) { aggr.ens_pd.clear(); - aggr.crps_climo_na.clear(); + aggr.crpscl_emp_na.clear(); + aggr.crpscl_gaus_na.clear(); aggr.hdr.clear(); m[key] = aggr; } @@ -2571,7 +2572,10 @@ void aggr_ecnt_lines(LineDataFile &f, STATAnalysisJob &job, // // Store the current statistics and weight (TOTAL column) // - m[key].ens_pd.crps_na.add(cur.crps); + m[key].ens_pd.crps_emp_na.add(cur.crps_emp); + m[key].crpscl_emp_na.add(cur.crpscl_emp); + m[key].ens_pd.crps_gaus_na.add(cur.crps_gaus); + m[key].crpscl_gaus_na.add(cur.crpscl_gaus); m[key].ens_pd.ign_na.add(cur.ign); m[key].ens_pd.var_na.add(square(cur.spread)); m[key].ens_pd.var_oerr_na.add(square(cur.spread_oerr)); @@ -2590,18 +2594,6 @@ void aggr_ecnt_lines(LineDataFile &f, STATAnalysisJob &job, bad_data_double : cur.rmse_oerr * cur.rmse_oerr)); - // - // Compute and store climatological CRPS - // - if(!is_bad_data(cur.crps) && !is_bad_data(cur.crpss) && - !is_eq(cur.crpss, 1.0)) { - crps_climo = cur.crps / (1.0 - cur.crpss); - m[key].crps_climo_na.add(crps_climo); - } - else { - m[key].crps_climo_na.add(bad_data_double); - } - // // Keep track of the unique header column entries // @@ -2627,17 +2619,20 @@ void aggr_ecnt_lines(LineDataFile &f, STATAnalysisJob &job, v = it->second.mse_oerr_na.wmean(it->second.ens_pd.wgt_na); it->second.ens_pd.rmse_oerr = (is_bad_data(v) ? bad_data_double : sqrt(v)); - crps_fcst = it->second.ens_pd.crps_na.wmean(it->second.ens_pd.wgt_na); - crps_climo = it->second.crps_climo_na.wmean(it->second.ens_pd.wgt_na); + crps_emp = it->second.ens_pd.crps_emp_na.wmean(it->second.ens_pd.wgt_na); + crpscl_emp = it->second.crpscl_emp_na.wmean(it->second.ens_pd.wgt_na); + crps_gaus = it->second.ens_pd.crps_gaus_na.wmean(it->second.ens_pd.wgt_na); + crpscl_gaus = it->second.crpscl_gaus_na.wmean(it->second.ens_pd.wgt_na); - // Compute aggregated CRPSS - if(!is_bad_data(crps_fcst) && !is_bad_data(crps_climo) && - !is_eq(crps_climo, 0.0)) { - it->second.ens_pd.crpss = (crps_climo - crps_fcst)/crps_climo; - } - else { - it->second.ens_pd.crpss = bad_data_double; - } + // Compute aggregated empirical CRPSS + it->second.ens_pd.crpss_emp = + (is_bad_data(crps_emp) || is_bad_data(crpscl_emp) && is_eq(crpscl_emp, 0.0) ? + bad_data_double : (crpscl_emp - crps_emp)/crpscl_emp); + + // Compute aggregated Gaussian CRPSS + it->second.ens_pd.crpss_gaus = + (is_bad_data(crps_gaus) || is_bad_data(crpscl_gaus) && is_eq(crpscl_gaus, 0.0) ? + bad_data_double : (crpscl_gaus - crps_gaus)/crpscl_gaus); } // end for it @@ -3008,7 +3003,7 @@ void aggr_orank_lines(LineDataFile &f, STATAnalysisJob &job, ORANKData cur; ConcatString key; int i, n_valid, n_bin; - double esum, esumsq, crps, ign, pit; + double esum, esumsq, crps_gaus, ign, pit; map::iterator it; // @@ -3102,9 +3097,10 @@ void aggr_orank_lines(LineDataFile &f, STATAnalysisJob &job, m[key].ens_pd.esumsq_na.add(esumsq); m[key].ens_pd.v_na.add(n_valid); - // Compute CRPS, IGN, and PIT for the current point - compute_crps_ign_pit(cur.obs, cur.ens_na, crps, ign, pit); - m[key].ens_pd.crps_na.add(crps); + // Compute Empirical and Gaussian CRPS, IGN, and PIT for the current point + m[key].ens_pd.crps_emp_na.add(compute_crps_emp(cur.obs, cur.ens_na)); + compute_crps_gaus_ign_pit(cur.obs, cur.ens_na, crps_gaus, ign, pit); + m[key].ens_pd.crps_gaus_na.add(crps_gaus); m[key].ens_pd.ign_na.add(ign); m[key].ens_pd.pit_na.add(pit); diff --git a/met/src/tools/core/stat_analysis/aggr_stat_line.h b/met/src/tools/core/stat_analysis/aggr_stat_line.h index 05f404cf6e..071645ade9 100644 --- a/met/src/tools/core/stat_analysis/aggr_stat_line.h +++ b/met/src/tools/core/stat_analysis/aggr_stat_line.h @@ -151,7 +151,7 @@ struct AggrISCInfo { struct AggrENSInfo { StatHdrInfo hdr; PairDataEnsemble ens_pd; - NumArray crps_climo_na; + NumArray crpscl_emp_na, crpscl_gaus_na; NumArray me_na, mse_na, me_oerr_na, mse_oerr_na; }; diff --git a/met/src/tools/core/stat_analysis/parse_stat_line.cc b/met/src/tools/core/stat_analysis/parse_stat_line.cc index e2d10201a0..365be05c90 100644 --- a/met/src/tools/core/stat_analysis/parse_stat_line.cc +++ b/met/src/tools/core/stat_analysis/parse_stat_line.cc @@ -347,10 +347,15 @@ void parse_ecnt_line(STATLine &l, ECNTData &e_data) { e_data.total = atoi(l.get_item("TOTAL")); e_data.n_ens = atof(l.get_item("N_ENS")); - e_data.crps = atof(l.get_item("CRPS")); - e_data.crpss = atof(l.get_item("CRPSS")); - e_data.ign = atof(l.get_item("IGN")); + e_data.crps_emp = atof(l.get_item("CRPS_EMP")); + e_data.crpscl_emp = atof(l.get_item("CRPSCL_EMP")); + e_data.crpss_emp = atof(l.get_item("CRPSS_EMP")); + + e_data.crps_gaus = atof(l.get_item("CRPS")); + e_data.crpscl_gaus = atof(l.get_item("CRPSCL")); + e_data.crpss_gaus = atof(l.get_item("CRPSS")); + e_data.ign = atof(l.get_item("IGN")); e_data.me = atof(l.get_item("ME")); e_data.rmse = atof(l.get_item("RMSE")); e_data.spread = atof(l.get_item("SPREAD")); diff --git a/met/src/tools/core/stat_analysis/parse_stat_line.h b/met/src/tools/core/stat_analysis/parse_stat_line.h index d3c6cbdf91..4f6dafa43b 100644 --- a/met/src/tools/core/stat_analysis/parse_stat_line.h +++ b/met/src/tools/core/stat_analysis/parse_stat_line.h @@ -62,8 +62,9 @@ struct MPRData { // Ensemble continuous statistics (ECNT) data structure struct ECNTData { int total, n_ens; - double crps, crpss, ign; - double me, rmse, spread; + double crps_emp, crpscl_emp, crpss_emp; + double crps_gaus, crpscl_gaus, crpss_gaus; + double ign, me, rmse, spread; double me_oerr, rmse_oerr, spread_oerr; double spread_plus_oerr; }; diff --git a/test/hdr/met_10_0.hdr b/test/hdr/met_10_0.hdr index 1c3004e680..8355e557a3 100644 --- a/test/hdr/met_10_0.hdr +++ b/test/hdr/met_10_0.hdr @@ -17,7 +17,7 @@ PJC : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_L PRC : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL N_THRESH _VAR_ PSTD : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL N_THRESH BASER BASER_NCL BASER_NCU RELIABILITY RESOLUTION UNCERTAINTY ROC_AUC BRIER BRIER_NCL BRIER_NCU BRIERCL BRIERCL_NCL BRIERCL_NCU BSS BSS_SMPL _VAR_ ECLV : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL BASE N_PTS _VAR_ -ECNT : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL N_ENS CRPS CRPSS IGN ME RMSE SPREAD ME_OERR RMSE_OERR SPREAD_OERR SPREAD_PLUS_OERR +ECNT : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL N_ENS CRPS CRPSS IGN ME RMSE SPREAD ME_OERR RMSE_OERR SPREAD_OERR SPREAD_PLUS_OERR CRPSCL CRPS_EMP CRPSCL_EMP CRPSS_EMP RPS : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL N_PROB RPS_REL RPS_RES RPS_UNC RPS RPSS RPSS_SMPL RPS_COMP RHIST : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL CRPS IGN N_RANK CRPSS SPREAD _VAR_ PHIST : VERSION MODEL DESC FCST_LEAD FCST_VALID_BEG FCST_VALID_END OBS_LEAD OBS_VALID_BEG OBS_VALID_END FCST_VAR FCST_UNITS FCST_LEV OBS_VAR OBS_UNITS OBS_LEV OBTYPE VX_MASK INTERP_MTHD INTERP_PNTS FCST_THRESH OBS_THRESH COV_THRESH ALPHA LINE_TYPE TOTAL BIN_SIZE N_BIN _VAR_