Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Bugfix 1638 develop climo_cdf #1639

Merged
merged 3 commits into from Jan 26, 2021
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
2 changes: 1 addition & 1 deletion met/src/basic/vx_util/data_plane_util.cc
Expand Up @@ -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
Expand Down
35 changes: 33 additions & 2 deletions met/src/libcode/vx_statistics/pair_base.cc
Expand Up @@ -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
Expand Down
8 changes: 4 additions & 4 deletions met/src/tools/core/grid_stat/grid_stat.cc
Expand Up @@ -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);
}

Expand Down Expand Up @@ -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();
Expand Down