diff --git a/met/data/config/GridStatConfig_default b/met/data/config/GridStatConfig_default index 89dc90e156..2ee10600a2 100644 --- a/met/data/config/GridStatConfig_default +++ b/met/data/config/GridStatConfig_default @@ -107,6 +107,7 @@ climo_cdf = { cdf_bins = 1; center_bins = FALSE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/met/data/config/PointStatConfig_default b/met/data/config/PointStatConfig_default index 558e419726..3967152fe6 100644 --- a/met/data/config/PointStatConfig_default +++ b/met/data/config/PointStatConfig_default @@ -149,6 +149,7 @@ climo_cdf = { cdf_bins = 1; center_bins = FALSE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/met/data/config/SeriesAnalysisConfig_default b/met/data/config/SeriesAnalysisConfig_default index c725cb2a5c..0249671fb6 100644 --- a/met/data/config/SeriesAnalysisConfig_default +++ b/met/data/config/SeriesAnalysisConfig_default @@ -85,6 +85,12 @@ climo_stdev = { file_name = []; } +climo_cdf = { + cdf_bins = 1; + center_bins = FALSE; + direct_prob = FALSE; +} + //////////////////////////////////////////////////////////////////////////////// // @@ -111,8 +117,9 @@ mask = { } // -// Number of grid points to be processed concurrently. Set smaller to use -// less memory but increase the number of passes through the data. +// Number of grid points to be processed concurrently. Set smaller to use less +// memory but increase the number of passes through the data. If set <= 0, all +// grid points are processed concurrently. // block_size = 1024; diff --git a/met/docs/Users_Guide/config_options.rst b/met/docs/Users_Guide/config_options.rst index b56cb65d75..9f568932bd 100644 --- a/met/docs/Users_Guide/config_options.rst +++ b/met/docs/Users_Guide/config_options.rst @@ -1307,8 +1307,8 @@ over the "climo_mean" setting and then updating the "file_name" entry. The "climo_cdf" dictionary specifies how the the climatological mean ("climo_mean") and standard deviation ("climo_stdev") data are used to evaluate model performance relative to where the observation value falls -within the climatological distribution. This dictionary consists of 3 -entries: +within the climatological distribution. This dictionary consists of the +following entries: (1) The "cdf_bins" entry defines the climatological bins either as an integer @@ -1320,6 +1320,8 @@ entries: (3) The "write_bins" entry may be set to TRUE or FALSE. +(4) The "direct_prob" entry may be set to TRUE or FALSE. + MET uses the climatological mean and standard deviation to construct a normal PDF at each observation location. The total area under the PDF is 1, and the climatological CDF value is computed as the area of the PDF to the left of @@ -1348,11 +1350,11 @@ an even number of bins can only be uncentered. For example: 5 centered bins (cdf_bins = 5; center_bins = TRUE;) yields: 0.0, 0.125, 0.375, 0.625, 0.875, 1.0 -When multiple climatological bins are used, statistics are computed -separately for each bin, and the average of the statistics across those bins -is written to the output. When "write_bins" is true, the statistics for each -bin are also written to the output. The bin number is appended to the -contents of the VX_MASK output column. +When multiple climatological bins are used for Point-Stat and Grid-Stat, +statistics are computed separately for each bin, and the average of the +statistics across those bins is written to the output. When "write_bins" +is true, the statistics for each bin are also written to the output. +The bin number is appended to the contents of the VX_MASK output column. Setting the number of bins to 1 effectively disables this logic by grouping all pairs into a single climatological bin. @@ -1363,6 +1365,7 @@ all pairs into a single climatological bin. cdf_bins = 11; or an array of floats center_bins = TRUE; or FALSE write_bins = FALSE; or TRUE + direct_prob = FALSE; or TRUE } .. _climato_data: @@ -1381,7 +1384,18 @@ directly to compute Brier Skill Score (BSS). When "climo_mean" and "climo_stdev" are both set to non-probability fields, the MET tools use the mean, standard deviation, and observation event threshold to derive a normal approximation of the climatological -probabilities. Those derived probability values are used to compute BSS. +probabilities. + +The "direct_prob" option controls the derivation logic. When "direct_prob" is +true, the climatological probability is computed directly from the +climatological distribution at each point as the area to the left of +the event threshold value. For greater-than or greater-than-or-equal-to +thresholds, 1.0 minus the area is used. When "direct_prob" is false, the +"cdf_bins" values are sampled from climatological distribution. The probability +is computed as the proportion of those samples which meet the threshold criteria. +In this way, the number of bins impacts the resolution of the climatological +probabilities. These derived probability values are used to compute the +climatological Brier Score and Brier Skill Score. .. _mask_missing_flag: @@ -3593,6 +3607,8 @@ Computation may be memory intensive, especially for large grids. The "block_size" entry sets the number of grid points to be processed concurrently (i.e. in one pass through a time series). Smaller values require less memory but increase the number of passes through the data. +If set less than or equal to 0, it is automatically reset to the number +of grid points, and they are all processed concurrently. .. code-block:: none diff --git a/met/docs/Users_Guide/series-analysis.rst b/met/docs/Users_Guide/series-analysis.rst index bc9095dbf1..640d948012 100644 --- a/met/docs/Users_Guide/series-analysis.rst +++ b/met/docs/Users_Guide/series-analysis.rst @@ -123,8 +123,7 @@ ____________________ block_size = 1024; -Number of grid points to be processed concurrently. Set smaller to use less memory but increase the number of passes through the data. The amount of memory the Series-Analysis tool consumes is determined by the size of the grid, the length of the series, and the block_size entry defined above. The larger this entry is set the faster the tool will run, subject to the amount of memory available on the machine. - +Number of grid points to be processed concurrently. Set smaller to use less memory but increase the number of passes through the data. The amount of memory the Series-Analysis tool consumes is determined by the size of the grid, the length of the series, and the block_size entry defined above. The larger this entry is set the faster the tool will run, subject to the amount of memory available on the machine. If set less than or equal to 0, it is automatically reset to the number of grid points, and they are all processed concurrently. ____________________ diff --git a/met/src/basic/vx_config/config_constants.h b/met/src/basic/vx_config/config_constants.h index 01f03b3840..156c55881b 100644 --- a/met/src/basic/vx_config/config_constants.h +++ b/met/src/basic/vx_config/config_constants.h @@ -300,10 +300,11 @@ struct RegridInfo { // struct ClimoCDFInfo { - bool flag; // Flag to turn on/off climo CDF logic - int n_bin; // Number of climo CDF cdf bins - ThreshArray cdf_ta; // Array of CDF thresholds - bool write_bins; // Flag for writing the individual bins + bool flag; // Flag to turn on/off climo CDF logic + int n_bin; // Number of climo CDF cdf bins + ThreshArray cdf_ta; // Array of CDF thresholds + bool write_bins; // Flag for writing the individual bins + bool direct_prob; // Flag for the direct computation of probs ClimoCDFInfo(); void clear(); @@ -665,6 +666,7 @@ static const char conf_key_climo_cdf[] = "climo_cdf"; static const char conf_key_cdf_bins[] = "cdf_bins"; static const char conf_key_center_bins[] = "center_bins"; static const char conf_key_write_bins[] = "write_bins"; +static const char conf_key_direct_prob[] = "direct_prob"; static const char conf_key_time_interp_method[] = "time_interp_method"; static const char conf_key_day_interval[] = "day_interval"; static const char conf_key_hour_interval[] = "hour_interval"; diff --git a/met/src/basic/vx_config/config_util.cc b/met/src/basic/vx_config/config_util.cc index 1218e2940d..3d77ff1a8f 100644 --- a/met/src/basic/vx_config/config_util.cc +++ b/met/src/basic/vx_config/config_util.cc @@ -1498,6 +1498,7 @@ void ClimoCDFInfo::clear() { n_bin = 0; cdf_ta.clear(); write_bins = false; + direct_prob = false; } /////////////////////////////////////////////////////////////////////////////// @@ -1589,9 +1590,13 @@ ClimoCDFInfo parse_conf_climo_cdf(Dictionary *dict) { center = cdf_dict->lookup_bool(conf_key_center_bins); // Conf: write_bins - // Used by Grid-Stat and Point-Stat but not by Ensemble-Stat + // Used by Grid-Stat and Point-Stat + // Not used by Ensemble-Stat or Series-Analysis info.write_bins = cdf_dict->lookup_bool(conf_key_write_bins, false, false); + // Conf: direct_prob + info.direct_prob = cdf_dict->lookup_bool(conf_key_direct_prob, false, false); + // Check that at least one value is provided if(bins.n() == 0) { mlog << Error << "\nparse_conf_climo_cdf() -> " diff --git a/met/src/libcode/vx_nc_obs/met_point_data.cc b/met/src/libcode/vx_nc_obs/met_point_data.cc index 6cc15a7e81..f223c4acff 100644 --- a/met/src/libcode/vx_nc_obs/met_point_data.cc +++ b/met/src/libcode/vx_nc_obs/met_point_data.cc @@ -59,12 +59,6 @@ void MetPointData::init_from_scratch() { } -//////////////////////////////////////////////////////////////////////// - -//void MetPointData::allocate() { -// obs_data.allocate(); -//} - //////////////////////////////////////////////////////////////////////// void MetPointData::clear() { @@ -178,18 +172,10 @@ MetPointDataPython::MetPointDataPython() { MetPointDataPython::MetPointDataPython(MetPointDataPython &d) { init_from_scratch(); - //obs_data = d.get_point_obs_data(); - //header_data = d.get_header_data(); - obs_data->assign(*d.get_point_obs_data()); + obs_data = new MetPointObsData(); + MetPointObsData *from_obs_data = d.get_point_obs_data(); + if (from_obs_data) obs_data->assign(*from_obs_data); header_data.assign(*d.get_header_data()); -cout << " DEBUG HS MetPointData(MetPointData &d) is called \n"; -cout << " DEBUG HS MetPointData(MetPointData &d) &header_data.lat_array=" << &(header_data.lat_array) << "\n"; -cout << " DEBUG HS MetPointData(MetPointData &d) header_data.lat_array.n()=" << header_data.lat_array.n() << "\n"; -cout << " DEBUG HS MetPointData(MetPointData &d) header_data.lon_array.n()=" << header_data.lon_array.n() << "\n"; -cout << " DEBUG HS MetPointData(MetPointData &d) header_data.elv_array.n()=" << header_data.elv_array.n() << "\n"; -cout << " DEBUG HS MetPointData(MetPointData &d) header_data.typ_idx_array.n()=" << header_data.typ_idx_array.n() << "\n"; -cout << " DEBUG HS MetPointData(MetPointData &d) header_data.sid_idx_array.n()=" << header_data.sid_idx_array.n() << "\n"; -cout << " DEBUG HS MetPointData(MetPointData &d) header_data.vld_idx_array.n()=" << header_data.vld_idx_array.n() << "\n"; } @@ -199,17 +185,6 @@ MetPointDataPython::~MetPointDataPython() { clear(); } -//////////////////////////////////////////////////////////////////////// - -//void MetPointDataPython::init_from_scratch() { -// nobs = 0; -// nhdr = 0; -// qty_length = 0; -// -// use_var_id = false; -// use_arr_vars = false; -//} - //////////////////////////////////////////////////////////////////////// @@ -434,4 +409,3 @@ void MetPointHeader::reset_counters() { } /////////////////////////////////////////////////////////////////////////////// - diff --git a/met/src/libcode/vx_pointdata_python/pointdata_python.cc b/met/src/libcode/vx_pointdata_python/pointdata_python.cc index 846966a6f1..316ec0dfc2 100644 --- a/met/src/libcode/vx_pointdata_python/pointdata_python.cc +++ b/met/src/libcode/vx_pointdata_python/pointdata_python.cc @@ -84,8 +84,6 @@ mlog << Error << "\nMetPythonPointDataFile::operator=(const MetPythonPointDataFi exit ( 1 ); -return ( * this ); - } diff --git a/met/src/libcode/vx_statistics/compute_stats.cc b/met/src/libcode/vx_statistics/compute_stats.cc index 4d6a2046c1..4938565694 100644 --- a/met/src/libcode/vx_statistics/compute_stats.cc +++ b/met/src/libcode/vx_statistics/compute_stats.cc @@ -762,7 +762,7 @@ void compute_pctinfo(const PairDataPoint &pd, bool pstd_flag, // Use input climatological probabilities or derive them if(cmn_flag) { if(cprob_in) climo_prob = *cprob_in; - else climo_prob = derive_climo_prob(pd.cdf_info, + else climo_prob = derive_climo_prob(pd.cdf_info_ptr, pd.cmn_na, pd.csd_na, pct_info.othresh); } diff --git a/met/src/libcode/vx_statistics/ens_stats.cc b/met/src/libcode/vx_statistics/ens_stats.cc index ac9805a928..e69b74f786 100644 --- a/met/src/libcode/vx_statistics/ens_stats.cc +++ b/met/src/libcode/vx_statistics/ens_stats.cc @@ -522,7 +522,7 @@ void RPSInfo::set(const PairDataEnsemble &pd) { climo_pct.zero_out(); // Derive climatological probabilities - if(cmn_flag) climo_prob = derive_climo_prob(pd.cdf_info, + if(cmn_flag) climo_prob = derive_climo_prob(pd.cdf_info_ptr, pd.cmn_na, pd.csd_na, fthresh[i]); diff --git a/met/src/libcode/vx_statistics/pair_base.cc b/met/src/libcode/vx_statistics/pair_base.cc index a8f2f2ff83..5a41f95d1b 100644 --- a/met/src/libcode/vx_statistics/pair_base.cc +++ b/met/src/libcode/vx_statistics/pair_base.cc @@ -64,9 +64,11 @@ void PairBase::clear() { IsPointVx = false; mask_name.clear(); - mask_area_ptr = (MaskPlane *) 0; // Not allocated - mask_sid_ptr = (StringArray *) 0; // Not allocated - mask_llpnt_ptr = (MaskLatLon *) 0; // Not allocated + mask_area_ptr = (MaskPlane *) 0; // Not allocated + mask_sid_ptr = (StringArray *) 0; // Not allocated + mask_llpnt_ptr = (MaskLatLon *) 0; // Not allocated + + cdf_info_ptr = (const ClimoCDFInfo *) 0; // Not allocated msg_typ.clear(); msg_typ_vals.clear(); @@ -75,8 +77,6 @@ void PairBase::clear() { interp_mthd = InterpMthd_None; interp_shape = GridTemplateFactory::GridTemplate_None; - cdf_info.clear(); - o_na.clear(); x_na.clear(); y_na.clear(); @@ -114,9 +114,11 @@ void PairBase::erase() { IsPointVx = false; mask_name.erase(); - mask_area_ptr = (MaskPlane *) 0; // Not allocated - mask_sid_ptr = (StringArray *) 0; // Not allocated - mask_llpnt_ptr = (MaskLatLon *) 0; // Not allocated + mask_area_ptr = (MaskPlane *) 0; // Not allocated + mask_sid_ptr = (StringArray *) 0; // Not allocated + mask_llpnt_ptr = (MaskLatLon *) 0; // Not allocated + + cdf_info_ptr = (const ClimoCDFInfo *) 0; // Not allocated msg_typ.clear(); msg_typ_vals.clear(); @@ -124,8 +126,6 @@ void PairBase::erase() { interp_mthd = InterpMthd_None; interp_shape = GridTemplateFactory::GridTemplate_None; - cdf_info.clear(); - o_na.erase(); x_na.erase(); y_na.erase(); @@ -218,6 +218,15 @@ void PairBase::set_mask_llpnt_ptr(MaskLatLon *llpnt_ptr) { //////////////////////////////////////////////////////////////////////// +void PairBase::set_climo_cdf_info_ptr(const ClimoCDFInfo *info_ptr) { + + cdf_info_ptr = info_ptr; + + return; +} + +//////////////////////////////////////////////////////////////////////// + void PairBase::set_msg_typ(const char *c) { msg_typ = c; @@ -272,15 +281,6 @@ void PairBase::set_interp_shape(GridTemplateFactory::GridTemplates shape) { //////////////////////////////////////////////////////////////////////// -void PairBase::set_climo_cdf_info(const ClimoCDFInfo &info) { - - cdf_info = info; - - return; -} - -//////////////////////////////////////////////////////////////////////// - void PairBase::set_fcst_ut(unixtime ut){ fcst_ut = ut; @@ -1027,22 +1027,25 @@ bool set_climo_flag(const NumArray &f_na, const NumArray &c_na) { //////////////////////////////////////////////////////////////////////// -void derive_climo_vals(const ClimoCDFInfo &cdf_info, +void derive_climo_vals(const ClimoCDFInfo *cdf_info_ptr, double m, double s, NumArray &climo_vals) { // Initialize climo_vals.erase(); - // cdf_info.cdf_ta starts with >=0.0 and ends with >=1.0. + // Check for no work to do + if(!cdf_info_ptr) return; + + // cdf_info_ptr->cdf_ta starts with >=0.0 and ends with >=1.0. // The number of bins is the number of thresholds minus 1. // Check for bad mean value - if(is_bad_data(m) || cdf_info.cdf_ta.n() < 2) { + if(is_bad_data(m) || cdf_info_ptr->cdf_ta.n() < 2) { return; } // Single climo bin - else if(cdf_info.cdf_ta.n() == 2) { + else if(cdf_info_ptr->cdf_ta.n() == 2) { climo_vals.add(m); } // Check for bad standard deviation value @@ -1053,9 +1056,9 @@ void derive_climo_vals(const ClimoCDFInfo &cdf_info, else { // Skip the first and last thresholds - for(int i=1; icdf_ta.n()-1; i++) { climo_vals.add( - normal_cdf_inv(cdf_info.cdf_ta[i].get_value(), m, s)); + normal_cdf_inv(cdf_info_ptr->cdf_ta[i].get_value(), m, s)); } } @@ -1064,7 +1067,7 @@ void derive_climo_vals(const ClimoCDFInfo &cdf_info, //////////////////////////////////////////////////////////////////////// -NumArray derive_climo_prob(const ClimoCDFInfo &cdf_info, +NumArray derive_climo_prob(const ClimoCDFInfo *cdf_info_ptr, const NumArray &mn_na, const NumArray &sd_na, const SingleThresh &othresh) { int i, n_mn, n_sd; @@ -1090,21 +1093,74 @@ NumArray derive_climo_prob(const ClimoCDFInfo &cdf_info, // threshold else if(n_mn > 0 && n_sd > 0) { - // The first (>=0.0) and last (>=1.0) climo thresholds are omitted - mlog << Debug(4) - << "Deriving climatological probabilities for threshold " - << othresh.get_str() << " by sampling " << cdf_info.cdf_ta.n()-2 - << " values from the normal climatological distribution.\n"; - - // Compute the probability by sampling from the climo distribution - // and deriving the event frequency - for(i=0; idirect_prob) { + + mlog << Debug(4) + << "Deriving climatological probabilities for threshold " + << othresh.get_str() << " directly from the normal " + << "climatological distribution.\n"; + + // Directly derive the climatological probability + for(i=0; icdf_ta.n() == 0) { + mlog << Error << "\nderive_climo_prob() -> " + << "No climatological probability thresholds defined!\n\n"; + exit(1); + } + + // The first (>=0.0) and last (>=1.0) climo thresholds are omitted + mlog << Debug(4) + << "Deriving climatological probabilities for threshold " + << othresh.get_str() << " by sampling " + << cdf_info_ptr->cdf_ta.n()-2 + << " values from the normal climatological distribution.\n"; + + // Compute the probability by sampling from the climo distribution + // and deriving the event frequency + for(i=0; icdf_ta.n() == 2) { mlog << Debug(3) << "Computing ensemble statistics relative to the " << "climatological mean.\n"; } - else if(cmn_flag && csd_flag && cdf_info.cdf_ta.n() > 2) { + else if(cmn_flag && csd_flag && cdf_info_ptr && cdf_info_ptr->cdf_ta.n() > 2) { mlog << Debug(3) << "Computing ensemble statistics relative to a " - << cdf_info.cdf_ta.n() - 2 + << cdf_info_ptr->cdf_ta.n() - 2 << "-member climatological ensemble.\n"; } else { @@ -468,7 +468,7 @@ void PairDataEnsemble::compute_pair_vals(const gsl_rng *rng_ptr) { } // Derive ensemble from climo mean and standard deviation - derive_climo_vals(cdf_info, cmn_na[i], csd_na[i], cur_clm); + derive_climo_vals(cdf_info_ptr, cmn_na[i], csd_na[i], cur_clm); // Store empirical CRPS stats crps_emp_na.add(compute_crps_emp(o_na[i], cur_ens)); @@ -785,7 +785,7 @@ PairDataEnsemble PairDataEnsemble::subset_pairs_obs_thresh(const SingleThresh &o pd.ssvar_bin_size = ssvar_bin_size; pd.obs_error_entry = obs_error_entry; pd.obs_error_flag = obs_error_flag; - pd.cdf_info = cdf_info; + pd.cdf_info_ptr = cdf_info_ptr; bool cmn_flag = set_climo_flag(o_na, cmn_na); bool csd_flag = set_climo_flag(o_na, csd_na); @@ -1306,12 +1306,12 @@ void VxPairDataEnsemble::set_ens_size(int n) { //////////////////////////////////////////////////////////////////////// -void VxPairDataEnsemble::set_climo_cdf_info(const ClimoCDFInfo &info) { +void VxPairDataEnsemble::set_climo_cdf_info_ptr(const ClimoCDFInfo *info) { for(int i=0; i " + << "Ensemble post-processing should be moved to the " + << "Gen-Ens-Prod tool, which replaces the logic of the " + << "\"ens\" dictionary. Support for the \"ens\" dictionary " + << "will be deprecated and removed." << "\n\n"; + } + // Parse the ensemble field information for(i=0,max_n_thresh=0; imsg_typ_sfc); diff --git a/met/src/tools/core/grid_stat/grid_stat.cc b/met/src/tools/core/grid_stat/grid_stat.cc index e7055ed596..8549cbf2e8 100644 --- a/met/src/tools/core/grid_stat/grid_stat.cc +++ b/met/src/tools/core/grid_stat/grid_stat.cc @@ -1847,7 +1847,7 @@ void get_mask_points(const GridStatVxOpt &vx_opt, pd.erase(); // Store the climo CDF info - pd.set_climo_cdf_info(vx_opt.cdf_info); + pd.set_climo_cdf_info_ptr(&vx_opt.cdf_info); // Apply the mask the data fields or fill with default values apply_mask(*fcst_ptr, mask_mp, pd.f_na); diff --git a/met/src/tools/core/point_stat/point_stat.cc b/met/src/tools/core/point_stat/point_stat.cc index 48ffd21c6e..11a88e0391 100644 --- a/met/src/tools/core/point_stat/point_stat.cc +++ b/met/src/tools/core/point_stat/point_stat.cc @@ -1736,7 +1736,7 @@ void do_hira_ens(int i_vx, const PairDataPoint *pd_ptr) { hira_pd.clear(); hira_pd.extend(pd_ptr->n_obs); hira_pd.set_ens_size(gt->size()); - hira_pd.set_climo_cdf_info(conf_info.vx_opt[i_vx].cdf_info); + hira_pd.set_climo_cdf_info_ptr(&conf_info.vx_opt[i_vx].cdf_info); f_ens.extend(gt->size()); // Process each observation point diff --git a/met/src/tools/core/point_stat/point_stat_conf_info.cc b/met/src/tools/core/point_stat/point_stat_conf_info.cc index 6117a0f7fc..8b90544b44 100644 --- a/met/src/tools/core/point_stat/point_stat_conf_info.cc +++ b/met/src/tools/core/point_stat/point_stat_conf_info.cc @@ -972,7 +972,7 @@ void PointStatVxOpt::set_vx_pd(PointStatConfInfo *conf_info) { vx_pd.set_mpr_thresh(mpr_sa, mpr_ta); // Store the climo CDF info - vx_pd.set_climo_cdf_info(cdf_info); + vx_pd.set_climo_cdf_info_ptr(&cdf_info); // Store the surface message type group cs = surface_msg_typ_group_str; diff --git a/met/src/tools/core/series_analysis/series_analysis.cc b/met/src/tools/core/series_analysis/series_analysis.cc index b89eeff0bb..b7f9aa01b9 100644 --- a/met/src/tools/core/series_analysis/series_analysis.cc +++ b/met/src/tools/core/series_analysis/series_analysis.cc @@ -30,6 +30,7 @@ // 010 12/11/19 Halley Gotway Reorganize logic to support the use // of python embedding. // 011 05/28/21 Halley Gotway Add MCTS HSS_EC output. +// 012 01/20/22 Halley Gotway MET #2003 Add PSTD BRIERCL output. // //////////////////////////////////////////////////////////////////////// @@ -339,6 +340,9 @@ void process_grid(const Grid &fcst_grid, const Grid &obs_grid) { // Process masking regions conf_info.process_masks(grid); + // Set the block size, if needed + if(is_bad_data(conf_info.block_size)) conf_info.block_size = nxy; + // Compute the number of reads required n_reads = nint(ceil((double) nxy / conf_info.block_size)); @@ -666,20 +670,10 @@ void process_scores() { // Number of points skipped due to valid data threshold int n_skip_zero = 0; int n_skip_pos = 0; - - // Allocate space to store the pairs for each grid point - pd_ptr = new PairDataPoint [conf_info.block_size]; - for(i=0; imagic_str() << ".\n"; @@ -1845,7 +1855,11 @@ void store_stat_pstd(int n, const ConcatString &col, else if(c == "BRIER") { v = pct_info.brier.v; } else if(c == "BRIER_NCL") { v = pct_info.brier.v_ncl[i]; } else if(c == "BRIER_NCU") { v = pct_info.brier.v_ncu[i]; } + else if(c == "BRIERCL") { v = pct_info.briercl.v; } + else if(c == "BRIERCL_NCL") { v = pct_info.briercl.v_ncl[i]; } + else if(c == "BRIERCL_NCU") { v = pct_info.briercl.v_ncu[i]; } else if(c == "BSS") { v = pct_info.bss; } + else if(c == "BSS_SMPL") { v = pct_info.bss_smpl; } else { mlog << Error << "\nstore_stat_pstd() -> " << "unsupported column name requested \"" << c diff --git a/met/src/tools/core/series_analysis/series_analysis_conf_info.cc b/met/src/tools/core/series_analysis/series_analysis_conf_info.cc index 0920787d0f..5fa6608650 100644 --- a/met/src/tools/core/series_analysis/series_analysis_conf_info.cc +++ b/met/src/tools/core/series_analysis/series_analysis_conf_info.cc @@ -68,6 +68,7 @@ void SeriesAnalysisConfInfo::clear() { fcnt_ta.clear(); ocnt_ta.clear(); cnt_logic = SetLogic_None; + cdf_info.clear(); ci_alpha.clear(); boot_interval = BootIntervalType_None; boot_rep_prop = bad_data_double; @@ -340,11 +341,11 @@ void SeriesAnalysisConfInfo::process_config(GrdFileType ftype, // Conf: block_size block_size = conf.lookup_int(conf_key_block_size); + // Reset invalid block_sizes to bad data so that they if(block_size <= 0.0) { - mlog << Error << "\nSeriesAnalysisConfInfo::process_config() -> " - << "The \"" << conf_key_block_size << "\" parameter (" - << block_size << ") must be greater than 0.\n\n"; - exit(1); + block_size = bad_data_int; + mlog << Debug(3) << "Automatically setting the \"" + << conf_key_block_size << "\" parameter to the size of the grid.\n"; } // Conf: vld_thresh @@ -414,6 +415,9 @@ void SeriesAnalysisConfInfo::process_config(GrdFileType ftype, } // end if continuous + // Conf: climo_cdf + cdf_info = parse_conf_climo_cdf(&conf); + // Conf: ci_alpha ci_alpha = parse_conf_ci_alpha(&conf); diff --git a/met/src/tools/core/series_analysis/series_analysis_conf_info.h b/met/src/tools/core/series_analysis/series_analysis_conf_info.h index 339de45ef4..f0cd829eed 100644 --- a/met/src/tools/core/series_analysis/series_analysis_conf_info.h +++ b/met/src/tools/core/series_analysis/series_analysis_conf_info.h @@ -57,6 +57,8 @@ class SeriesAnalysisConfInfo { ThreshArray ocnt_ta; // Continuous obs thresholds SetLogic cnt_logic; // Continuous threshold field logic + ClimoCDFInfo cdf_info; // Climo CDF info + NumArray ci_alpha; // Alpha value for confidence intervals BootIntervalType boot_interval; // Bootstrap CI type double boot_rep_prop; // Bootstrap replicate proportion 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 3e66de254e..edee31714d 100644 --- a/met/src/tools/core/stat_analysis/aggr_stat_line.cc +++ b/met/src/tools/core/stat_analysis/aggr_stat_line.cc @@ -2204,9 +2204,6 @@ void aggr_mpr_lines(LineDataFile &f, STATAnalysisJob &job, // if(m.count(key) == 0) { - bool center = false; - aggr.pd.cdf_info.set_cdf_ta(nint(1.0/job.out_bin_size), center); - aggr.pd.f_na.clear(); aggr.pd.o_na.clear(); aggr.pd.cmn_na.clear(); @@ -2226,8 +2223,14 @@ void aggr_mpr_lines(LineDataFile &f, STATAnalysisJob &job, aggr.obs_var = cur.obs_var; aggr.hdr.clear(); + bool center = false; + aggr.cdf_info.set_cdf_ta(nint(1.0/job.out_bin_size), center); + m[key] = aggr; + // Set the pointer after storing in the map + m[key].pd.set_climo_cdf_info_ptr(&m[key].cdf_info); + mlog << Debug(3) << "[Case " << m.size() << "] Added new case for key \"" << key << "\".\n"; @@ -3084,9 +3087,9 @@ void aggr_orank_lines(LineDataFile &f, STATAnalysisJob &job, // Add a new map entry, if necessary // if(m.count(key) == 0) { + aggr.clear(); - bool center = false; - aggr.ens_pd.cdf_info.set_cdf_ta(nint(1.0/job.out_bin_size), center); + aggr.ens_pd.obs_error_flag = !is_bad_data(cur.ens_mean_oerr); aggr.ens_pd.set_ens_size(cur.n_ens); aggr.ens_pd.extend(cur.total); @@ -3095,7 +3098,15 @@ void aggr_orank_lines(LineDataFile &f, STATAnalysisJob &job, n_bin = ceil(1.0/aggr.ens_pd.phist_bin_size); for(i=0; i= debug_level) { @@ -1995,11 +1996,16 @@ void process_pbfile(int i_pb) { } if(i_msg <= 0) { - mlog << Warning << "\n" << method_name - << ((n_derived_obs > 0) ? "Saved the derived variables only. " : " ") - << "No " << (is_prepbufr ? "PrepBufr" : "Bufr") - << " messages retained from file: " - << pbfile[i_pb] << "\n\n"; + if (n_derived_obs > 0) + mlog << Debug(3) << method_name + << "Saved the derived variables only. No " << (is_prepbufr ? "PrepBufr" : "Bufr") + << " messages retained from file: " + << pbfile[i_pb] << "\n"; + else + mlog << Warning << "\n" << method_name + << "No " << (is_prepbufr ? "PrepBufr" : "Bufr") + << " messages retained from file: " + << pbfile[i_pb] << "\n\n"; } return; diff --git a/test/config/GridStatConfig_climo_WMO b/test/config/GridStatConfig_climo_WMO index e23f74a2be..f16d03b47e 100644 --- a/test/config/GridStatConfig_climo_WMO +++ b/test/config/GridStatConfig_climo_WMO @@ -130,6 +130,7 @@ climo_cdf = { cdf_bins = 3; center_bins = TRUE; write_bins = FALSE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/GridStatConfig_climo_prob b/test/config/GridStatConfig_climo_prob index a97291af88..085a8135d9 100644 --- a/test/config/GridStatConfig_climo_prob +++ b/test/config/GridStatConfig_climo_prob @@ -140,6 +140,7 @@ climo_cdf = { cdf_bins = [ 0.0, 0.4, 0.6, 1.0 ]; center_bins = FALSE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/GridStatConfig_mpr_thresh b/test/config/GridStatConfig_mpr_thresh index 268e24cbc6..b2930f027b 100644 --- a/test/config/GridStatConfig_mpr_thresh +++ b/test/config/GridStatConfig_mpr_thresh @@ -128,6 +128,7 @@ climo_cdf = { cdf_bins = 1; center_bins = FALSE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/PointStatConfig_airnow b/test/config/PointStatConfig_airnow index 80df21807a..5c424a3cb6 100644 --- a/test/config/PointStatConfig_airnow +++ b/test/config/PointStatConfig_airnow @@ -126,6 +126,7 @@ climo_cdf = { cdf_bins = 1; center_bins = FALSE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/PointStatConfig_climo b/test/config/PointStatConfig_climo index 041642d1a6..207d251ff5 100644 --- a/test/config/PointStatConfig_climo +++ b/test/config/PointStatConfig_climo @@ -173,6 +173,7 @@ climo_cdf = { cdf_bins = 1; center_bins = FALSE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/PointStatConfig_climo_WMO b/test/config/PointStatConfig_climo_WMO index 599797cde3..c801f581ac 100644 --- a/test/config/PointStatConfig_climo_WMO +++ b/test/config/PointStatConfig_climo_WMO @@ -121,6 +121,7 @@ climo_cdf = { cdf_bins = 3; center_bins = TRUE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/PointStatConfig_climo_prob b/test/config/PointStatConfig_climo_prob index e57f5e003f..ca4cadcbc9 100644 --- a/test/config/PointStatConfig_climo_prob +++ b/test/config/PointStatConfig_climo_prob @@ -123,6 +123,7 @@ climo_cdf = { cdf_bins = 11; center_bins = TRUE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/PointStatConfig_mpr_thresh b/test/config/PointStatConfig_mpr_thresh index 6da3ee7c63..f1559f1d3f 100644 --- a/test/config/PointStatConfig_mpr_thresh +++ b/test/config/PointStatConfig_mpr_thresh @@ -115,6 +115,7 @@ climo_cdf = { cdf_bins = 1; center_bins = FALSE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/PointStatConfig_python b/test/config/PointStatConfig_python index 2e97f43f29..8d70867143 100644 --- a/test/config/PointStatConfig_python +++ b/test/config/PointStatConfig_python @@ -112,6 +112,7 @@ climo_cdf = { cdf_bins = 1; center_bins = FALSE; write_bins = TRUE; + direct_prob = FALSE; } //////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/SeriesAnalysisConfig b/test/config/SeriesAnalysisConfig index 8c5e965230..9d28ce3925 100644 --- a/test/config/SeriesAnalysisConfig +++ b/test/config/SeriesAnalysisConfig @@ -53,6 +53,39 @@ obs = { //////////////////////////////////////////////////////////////////////////////// +// +// Climatology data +// +climo_mean = { + + file_name = []; + field = []; + + regrid = { + method = NEAREST; + width = 1; + vld_thresh = 0.5; + shape = SQUARE; + } + + time_interp_method = DW_MEAN; + day_interval = 31; + hour_interval = 6; +} + +climo_stdev = climo_mean; +climo_stdev = { + file_name = []; +} + +climo_cdf = { + cdf_bins = 1; + center_bins = FALSE; + direct_prob = FALSE; +} + +//////////////////////////////////////////////////////////////////////////////// + // // Confidence interval settings // @@ -77,10 +110,11 @@ mask = { } // -// Number of grid points to be processed concurrently. Set smaller to use -// less memory but increase the number of passes through the data. +// Number of grid points to be processed concurrently. Set smaller to use less +// memory but increase the number of passes through the data. If set <= 0, all +// grid points are processed concurrently. // -block_size = 10000; +block_size = 169*154; // // Ratio of valid matched pairs to compute statistics for a grid point diff --git a/test/config/SeriesAnalysisConfig_climo b/test/config/SeriesAnalysisConfig_climo index a1ceb41577..78d9e7449f 100644 --- a/test/config/SeriesAnalysisConfig_climo +++ b/test/config/SeriesAnalysisConfig_climo @@ -82,6 +82,12 @@ climo_stdev = { file_name = [ ${CLIMO_STDEV_FILE_LIST} ]; } +climo_cdf = { + cdf_bins = 1; + center_bins = FALSE; + direct_prob = FALSE; +} + //////////////////////////////////////////////////////////////////////////////// // @@ -108,10 +114,11 @@ mask = { } // -// Number of grid points to be processed concurrently. Set smaller to use -// less memory but increase the number of passes through the data. +// Number of grid points to be processed concurrently. Set smaller to use less +// memory but increase the number of passes through the data. If set <= 0, all +// grid points are processed concurrently. // -block_size = 150000; +block_size = 0; // // Ratio of valid matched pairs to compute statistics for a grid point diff --git a/test/config/SeriesAnalysisConfig_climo_prob b/test/config/SeriesAnalysisConfig_climo_prob new file mode 100644 index 0000000000..ee023bf6c4 --- /dev/null +++ b/test/config/SeriesAnalysisConfig_climo_prob @@ -0,0 +1,164 @@ +//////////////////////////////////////////////////////////////////////////////// +// +// Series-Analysis configuration file. +// +// For additional information, please see the MET User's Guide. +// +//////////////////////////////////////////////////////////////////////////////// + +// +// Output model name to be written +// +model = "SREF"; + +// +// Output description to be written +// +desc = "NA"; + +// +// Output observation type to be written +// +obtype = "GFSANL"; + +//////////////////////////////////////////////////////////////////////////////// + +// +// Verification grid +// +regrid = { + to_grid = FCST; + method = BILIN; + width = 2; + vld_thresh = 0.5; + shape = SQUARE; +} + +//////////////////////////////////////////////////////////////////////////////// + +censor_thresh = []; +censor_val = []; +cat_thresh = []; +cnt_thresh = [ NA ]; +cnt_logic = UNION; + +// +// Forecast and observation fields to be verified +// +fcst = { + cat_thresh = ==0.1; + + field = [ + { + name = "PROB"; level = "Z2"; + prob = { name = "TMP"; thresh_hi = 273.0; }; + } + ]; +} +obs = { + cat_thresh = <273.0; + + field = [ + { + name = "TMP"; level = "Z2"; + } + ]; +} + +//////////////////////////////////////////////////////////////////////////////// + +// +// Climatology data +// +climo_mean = obs; +climo_mean = { + + file_name = [ ${CLIMO_MEAN_FILE_LIST} ]; + + regrid = { + method = BILIN; + width = 2; + vld_thresh = 0.5; + } + + time_interp_method = DW_MEAN; + day_interval = ${DAY_INTERVAL}; + hour_interval = ${HOUR_INTERVAL}; +} + +climo_stdev = climo_mean; +climo_stdev = { + file_name = [ ${CLIMO_STDEV_FILE_LIST} ]; +} + +climo_cdf = { + cdf_bins = 1; + center_bins = FALSE; + direct_prob = TRUE; +} + +//////////////////////////////////////////////////////////////////////////////// + +// +// Confidence interval settings +// +ci_alpha = [ 0.05 ]; + +boot = { + interval = PCTILE; + rep_prop = 1.0; + n_rep = 0; + rng = "mt19937"; + seed = "1"; +} + +//////////////////////////////////////////////////////////////////////////////// + +// +// Verification masking regions +// +mask = { + grid = ""; + poly = ""; +} + +// +// Number of grid points to be processed concurrently. Set smaller to use less +// memory but increase the number of passes through the data. If set <= 0, all +// grid points are processed concurrently. +// +block_size = -9999; + +// +// Ratio of valid matched pairs to compute statistics for a grid point +// +vld_thresh = 0.5; + +//////////////////////////////////////////////////////////////////////////////// + +// +// Statistical output types +// +output_stats = { + fho = [ ]; + ctc = [ ]; + cts = [ ]; + mctc = [ ]; + mcts = [ ]; + cnt = [ ]; + sl1l2 = [ ]; + sal1l2 = [ ]; + pct = [ ]; + pstd = [ "TOTAL", "ROC_AUC", "BRIER", "BRIERCL", "BSS", "BSS_SMPL" ]; + pjc = [ ]; + prc = [ ]; +} + +//////////////////////////////////////////////////////////////////////////////// + +hss_ec_value = NA; +rank_corr_flag = FALSE; +tmp_dir = "/tmp"; +version = "V10.1.0"; + +//////////////////////////////////////////////////////////////////////////////// diff --git a/test/config/SeriesAnalysisConfig_conditional b/test/config/SeriesAnalysisConfig_conditional index a4b1447a1c..b4becd6656 100644 --- a/test/config/SeriesAnalysisConfig_conditional +++ b/test/config/SeriesAnalysisConfig_conditional @@ -77,8 +77,9 @@ mask = { } // -// Number of grid points to be processed concurrently. Set smaller to use -// less memory but increase the number of passes through the data. +// Number of grid points to be processed concurrently. Set smaller to use less +// memory but increase the number of passes through the data. If set <= 0, all +// grid points are processed concurrently. // block_size = 10000; diff --git a/test/config/SeriesAnalysisConfig_python b/test/config/SeriesAnalysisConfig_python index 945a1a56b6..a33b237ab3 100644 --- a/test/config/SeriesAnalysisConfig_python +++ b/test/config/SeriesAnalysisConfig_python @@ -82,6 +82,12 @@ climo_stdev = { file_name = []; } +climo_cdf = { + cdf_bins = 1; + center_bins = FALSE; + direct_prob = FALSE; +} + //////////////////////////////////////////////////////////////////////////////// // @@ -108,8 +114,9 @@ mask = { } // -// Number of grid points to be processed concurrently. Set smaller to use -// less memory but increase the number of passes through the data. +// Number of grid points to be processed concurrently. Set smaller to use less +// memory but increase the number of passes through the data. If set <= 0, all +// grid points are processed concurrently. // block_size = 24000; diff --git a/test/xml/unit_climatology_1.0deg.xml b/test/xml/unit_climatology_1.0deg.xml index 8fad5f3be1..6f458f13c8 100644 --- a/test/xml/unit_climatology_1.0deg.xml +++ b/test/xml/unit_climatology_1.0deg.xml @@ -190,6 +190,55 @@ + + echo "&DATA_DIR_MODEL;/grib2/sref_pr/sref_prob_2012040821_F003.grib2 \ + &DATA_DIR_MODEL;/grib2/sref_pr/sref_prob_2012040821_F009.grib2 \ + &DATA_DIR_MODEL;/grib2/sref_pr/sref_prob_2012040821_F015.grib2 \ + &DATA_DIR_MODEL;/grib2/sref_pr/sref_prob_2012040821_F021.grib2 \ + &DATA_DIR_MODEL;/grib2/sref_pr/sref_prob_2012040821_F027.grib2 \ + &DATA_DIR_MODEL;/grib2/sref_pr/sref_prob_2012040821_F033.grib2 \ + &DATA_DIR_MODEL;/grib2/sref_pr/sref_prob_2012040821_F039.grib2 \ + &DATA_DIR_MODEL;/grib2/sref_pr/sref_prob_2012040821_F045.grib2" \ + > &OUTPUT_DIR;/climatology_1.0deg/input_fcst_file_list; \ + echo "&DATA_DIR_MODEL;/grib2/gfsanl/gfsanl_4_20120409_0000_000.grb2 \ + &DATA_DIR_MODEL;/grib2/gfsanl/gfsanl_4_20120409_0600_000.grb2 \ + &DATA_DIR_MODEL;/grib2/gfsanl/gfsanl_4_20120409_1200_000.grb2 \ + &DATA_DIR_MODEL;/grib2/gfsanl/gfsanl_4_20120409_1800_000.grb2 \ + &DATA_DIR_MODEL;/grib2/gfsanl/gfsanl_4_20120410_0000_000.grb2 \ + &DATA_DIR_MODEL;/grib2/gfsanl/gfsanl_4_20120410_0600_000.grb2 \ + &DATA_DIR_MODEL;/grib2/gfsanl/gfsanl_4_20120410_1200_000.grb2 \ + &DATA_DIR_MODEL;/grib2/gfsanl/gfsanl_4_20120410_1800_000.grb2" \ + > &OUTPUT_DIR;/climatology_1.0deg/input_obs_file_list; \ + &MET_BIN;/series_analysis + + DAY_INTERVAL 1 + HOUR_INTERVAL 6 + CLIMO_MEAN_FILE_LIST + "&DATA_DIR_CLIMO;/NCEP_NCAR_40YR_1.0deg/cmean_1d.19590409", + "&DATA_DIR_CLIMO;/NCEP_NCAR_40YR_1.0deg/cmean_1d.19590410", + "&DATA_DIR_CLIMO;/NCEP_NCAR_40YR_1.0deg/cmean_1d.19590411" + + + CLIMO_STDEV_FILE_LIST + "&DATA_DIR_CLIMO;/NCEP_NCAR_40YR_1.0deg/cstdv_1d.19590409", + "&DATA_DIR_CLIMO;/NCEP_NCAR_40YR_1.0deg/cstdv_1d.19590410", + "&DATA_DIR_CLIMO;/NCEP_NCAR_40YR_1.0deg/cstdv_1d.19590411" + + + + \ + -fcst &OUTPUT_DIR;/climatology_1.0deg/input_fcst_file_list \ + -obs &OUTPUT_DIR;/climatology_1.0deg/input_obs_file_list \ + -paired \ + -out &OUTPUT_DIR;/climatology_1.0deg/series_analysis_PROB_CLIMO_1.0DEG.nc \ + -config &CONFIG_DIR;/SeriesAnalysisConfig_climo_prob \ + -v 2 + + + &OUTPUT_DIR;/climatology_1.0deg/series_analysis_PROB_CLIMO_1.0DEG.nc + + + echo "&DATA_DIR_MODEL;/grib1/arw-fer-gep1/arw-fer-gep1_2012040912_F024.grib \ &DATA_DIR_MODEL;/grib1/arw-fer-gep5/arw-fer-gep5_2012040912_F024.grib \