diff --git a/met/data/config/TCGenConfig_default b/met/data/config/TCGenConfig_default index 0de6a27844..839d8bd21d 100644 --- a/met/data/config/TCGenConfig_default +++ b/met/data/config/TCGenConfig_default @@ -213,9 +213,9 @@ ops_method_flag = TRUE; //////////////////////////////////////////////////////////////////////////////// // -// Genesis probability thresholds +// Probability of genesis thresholds // -genesis_prob_thresh = ==0.25; +prob_genesis_thresh = ==0.25; // // Confidence interval alpha value diff --git a/met/src/basic/vx_config/config_constants.h b/met/src/basic/vx_config/config_constants.h index 3082eb69d7..01c598bd72 100644 --- a/met/src/basic/vx_config/config_constants.h +++ b/met/src/basic/vx_config/config_constants.h @@ -1097,7 +1097,7 @@ static const char conf_key_ops_hit_window[] = "ops_hit_window"; static const char conf_key_discard_init_post_genesis_flag[] = "discard_init_post_genesis_flag"; static const char conf_key_dev_method_flag[] = "dev_method_flag"; static const char conf_key_ops_method_flag[] = "ops_method_flag"; -static const char conf_key_genesis_prob_thresh[] = "genesis_prob_thresh"; +static const char conf_key_prob_genesis_thresh[] = "prob_genesis_thresh"; static const char conf_key_fcst_fy_oy[] = "fcst_fy_oy"; static const char conf_key_fcst_fy_on[] = "fcst_fy_on"; static const char conf_key_fcst_tracks[] = "fcst_tracks"; diff --git a/met/src/tools/tc_utils/tc_gen/tc_gen.cc b/met/src/tools/tc_utils/tc_gen/tc_gen.cc index f2a26967d8..dcbd87f7c6 100644 --- a/met/src/tools/tc_utils/tc_gen/tc_gen.cc +++ b/met/src/tools/tc_utils/tc_gen/tc_gen.cc @@ -85,6 +85,11 @@ static void get_genesis_pairs (const TCGenVxOpt &, static void do_genesis_ctc (const TCGenVxOpt &, PairDataGenesis &, GenCTCInfo &); +static void do_probgen_pct (const TCGenVxOpt &, + const ProbInfoArray &, + const GenesisInfoArray &, + const TrackInfoArray &, + ProbGenPCTInfo &); static int find_genesis_match (const GenesisInfo &, const GenesisInfoArray &, @@ -97,7 +102,8 @@ static void setup_nc_file (); static void write_ctc_stats (const PairDataGenesis &, GenCTCInfo &); -static void write_pct_stats (int, PCTInfo &, PCTInfo &); +static void write_pct_stats (ProbGenPCTInfo &); + static void write_genmpr_row (StatHdrColumns &, const PairDataGenesis &, STATOutputType, @@ -271,7 +277,7 @@ void score_track_genesis(const GenesisInfoArray &best_ga, map model_ga_map; map::iterator it; PairDataGenesis pairs; - GenCTCInfo ctc_info; + GenCTCInfo genesis_ctc; // Get the list of genesis track files get_atcf_files(genesis_source, genesis_model_suffix, atcf_gen_reg_exp, @@ -323,14 +329,14 @@ void score_track_genesis(const GenesisInfoArray &best_ga, } // end for j - // Process the genesis events for each model. + // Process the genesis events for each model for(j=0,it=model_ga_map.begin(); it!=model_ga_map.end(); it++,j++) { // Initialize - ctc_info.clear(); - ctc_info.Model = it->first; - ctc_info.set_vx_opt(&conf_info.VxOpt[i], - &conf_info.NcOutGrid); + genesis_ctc.clear(); + genesis_ctc.Model = it->first; + genesis_ctc.set_vx_opt(&conf_info.VxOpt[i], + &conf_info.NcOutGrid); mlog << Debug(2) << "[Filter " << i+1 << " (" << conf_info.VxOpt[i].Desc @@ -346,14 +352,14 @@ void score_track_genesis(const GenesisInfoArray &best_ga, best_ga, oper_ta, pairs); // Do the categorical verification - do_genesis_ctc(conf_info.VxOpt[i], pairs, ctc_info); + do_genesis_ctc(conf_info.VxOpt[i], pairs, genesis_ctc); // Write the statistics output - write_ctc_stats(pairs, ctc_info); + write_ctc_stats(pairs, genesis_ctc); // Write NetCDF output fields if(!conf_info.VxOpt[i].NcInfo.all_false()) { - write_nc(ctc_info); + write_nc(genesis_ctc); } } // end for j @@ -367,16 +373,13 @@ void score_track_genesis(const GenesisInfoArray &best_ga, void score_genesis_prob(const GenesisInfoArray &best_ga, const TrackInfoArray &oper_ta) { - int i, j, k, l, prob_time; - double prob_value; - IntArray prob_times; + int i, j; StringArray edeck_files, edeck_files_model_suffix; ProbInfoArray fcst_pa, empty_pa; - ConcatString model, cs; + ConcatString model; map model_pa_map; - map::iterator model_it; - map pct_dev_map, pct_ops_map; - PCTInfo empty_pct; + map::iterator it; + ProbGenPCTInfo probgen_pct; // Get the list of EDECK files get_atcf_files(edeck_source, edeck_model_suffix, atcf_reg_exp, @@ -397,8 +400,8 @@ void score_genesis_prob(const GenesisInfoArray &best_ga, // Initialize model_pa_map.clear(); - // Subset the forecast genesis probabilities - // Organize them into a map by model name + // Organize the forecast genesis probabilities + // into a map by model name for(j=0; jfirst - << " model, comparing " << model_it->second.n_prob_gen() - << " genesis forecasts to " << best_ga.n() << " " + << ") " << ": Model " << j+1 << "] " << "For " << it->first + << " model, comparing " << it->second.n_prob_gen() + << " probability of genesis forecasts to " << best_ga.n() << " " << conf_info.BestEventInfo.Technique << " and " << oper_ta.n() << " " << conf_info.OperTechnique << " tracks.\n"; - // Loop over the probability of genesis events - for(k=0; ksecond.n_prob_gen(); k++) { - - // JHG search for a match! - // For now, just assume it's a hit. - // Also need to define SHC values. - bool dev_match = true; - bool ops_match = true; - - // Loop over the individual probabilities - for(l=0; lsecond.prob_gen(k).n_prob(); l++) { - - // Current lead time - prob_time = nint(model_it->second.prob_gen(k).prob_item(l)); - prob_value = model_it->second.prob_gen(k).prob(l) / 100.0; - - // Add a new map entries for this lead time, if necessary - if(!prob_times.has(prob_time)) { - prob_times.add(prob_time); - pct_dev_map[prob_time] = empty_pct; - pct_ops_map[prob_time] = empty_pct; - } - - // Increment counts - if(dev_match) pct_dev_map[prob_time].pct.inc_event (prob_value); - else pct_dev_map[prob_time].pct.inc_nonevent(prob_value); - if(ops_match) pct_ops_map[prob_time].pct.inc_event (prob_value); - else pct_ops_map[prob_time].pct.inc_nonevent(prob_value); - - } // end for l - } // end for k - - // Write probabilistic output - for(k=0; ksecond, + best_ga, oper_ta, probgen_pct); + + // Write the statistics output + write_pct_stats(probgen_pct); } // end for j - } // end for i n_vx + } // end for i return; } @@ -758,6 +718,51 @@ void do_genesis_ctc(const TCGenVxOpt &vx_opt, return; } +//////////////////////////////////////////////////////////////////////// + +void do_probgen_pct(const TCGenVxOpt &vx_opt, + const ProbInfoArray &model_pa, + const GenesisInfoArray &best_ga, + const TrackInfoArray &oper_ta, + ProbGenPCTInfo &pgi) { + int i, j, prob_lead; + double prob_value; + + // Initialize + pgi.clear(); + pgi.set_vx_opt(&vx_opt); + + // Score each of the probability forecasts + for(i=0; iDesc.c_str()); + shc.set_obtype(conf_info.BestEventInfo.Technique.c_str()); + shc.set_mask(pgi.VxOpt->VxMaskName.empty() ? + na_str : pgi.VxOpt->VxMaskName.c_str()); - if(conf_info.VxOpt[i_vx].DevFlag) { - write_pct_row(shc, pct_dev, - conf_info.VxOpt[i_vx].OutputMap[txt_file_type[i_pct]], - 1, 1, stat_at, i_stat_row, - txt_at[i_pct], i_txt_row[i_pct]); - } - if(conf_info.VxOpt[i_vx].OpsFlag) { - write_pct_row(shc, pct_ops, - conf_info.VxOpt[i_vx].OutputMap[txt_file_type[i_pct]], - 1, 1, stat_at, i_stat_row, - txt_at[i_pct], i_txt_row[i_pct]); - } - } + // Write results for each lead time + for(i=0; ioutput_map(stat_pct) != STATOutputType_None) { + + if(pgi.VxOpt->DevFlag) { + shc.set_fcst_var(probgen_dev_name); + shc.set_obs_var (probgen_dev_name); + write_pct_row(shc, pgi.PCTDev[lead_hr], + pgi.VxOpt->output_map(stat_pct), + 1, 1, stat_at, i_stat_row, + txt_at[i_pct], i_txt_row[i_pct]); + } + if(pgi.VxOpt->OpsFlag) { + shc.set_fcst_var(probgen_ops_name); + shc.set_obs_var (probgen_ops_name); + write_pct_row(shc, pgi.PCTOps[lead_hr], + pgi.VxOpt->output_map(stat_pct), + 1, 1, stat_at, i_stat_row, + txt_at[i_pct], i_txt_row[i_pct]); + } } - if(conf_info.VxOpt[i_vx].OpsFlag) { - pct_ops.compute_stats(); - pct_ops.compute_ci(); - write_pstd_row(shc, pct_ops, - conf_info.VxOpt[i_vx].OutputMap[txt_file_type[i_pstd]], - 1, 1, stat_at, i_stat_row, - txt_at[i_pstd], i_txt_row[i_pstd]); + + // Write PSTD output + if(pgi.VxOpt->output_map(stat_pstd) != STATOutputType_None) { + + if(pgi.VxOpt->DevFlag) { + pgi.PCTDev[lead_hr].compute_stats(); + pgi.PCTDev[lead_hr].compute_ci(); + shc.set_fcst_var(probgen_dev_name); + shc.set_obs_var (probgen_dev_name); + write_pstd_row(shc, pgi.PCTDev[lead_hr], + pgi.VxOpt->output_map(stat_pstd), + 1, 1, stat_at, i_stat_row, + txt_at[i_pstd], i_txt_row[i_pstd]); + } + if(pgi.VxOpt->OpsFlag) { + pgi.PCTOps[lead_hr].compute_stats(); + pgi.PCTOps[lead_hr].compute_ci(); + shc.set_fcst_var(probgen_ops_name); + shc.set_obs_var (probgen_ops_name); + write_pstd_row(shc, pgi.PCTOps[lead_hr], + pgi.VxOpt->output_map(stat_pstd) , + 1, 1, stat_at, i_stat_row, + txt_at[i_pstd], i_txt_row[i_pstd]); + } } - } - // Write PJC output - if(conf_info.VxOpt[i_vx].OutputMap[txt_file_type[i_pjc]] != STATOutputType_None) { + // Write PJC output + if(pgi.VxOpt->output_map(stat_pjc) != STATOutputType_None) { - if(conf_info.VxOpt[i_vx].DevFlag) { - write_pct_row(shc, pct_dev, - conf_info.VxOpt[i_vx].OutputMap[txt_file_type[i_pjc]], - 1, 1, stat_at, i_stat_row, - txt_at[i_pjc], i_txt_row[i_pjc]); - } - if(conf_info.VxOpt[i_vx].OpsFlag) { - write_pct_row(shc, pct_ops, - conf_info.VxOpt[i_vx].OutputMap[txt_file_type[i_pjc]], - 1, 1, stat_at, i_stat_row, - txt_at[i_pjc], i_txt_row[i_pjc]); + if(pgi.VxOpt->DevFlag) { + shc.set_fcst_var(probgen_dev_name); + shc.set_obs_var (probgen_dev_name); + write_pct_row(shc, pgi.PCTDev[lead_hr], + pgi.VxOpt->output_map(stat_pjc), + 1, 1, stat_at, i_stat_row, + txt_at[i_pjc], i_txt_row[i_pjc]); + } + if(pgi.VxOpt->OpsFlag) { + shc.set_fcst_var(probgen_ops_name); + shc.set_obs_var (probgen_ops_name); + write_pct_row(shc, pgi.PCTOps[lead_hr], + pgi.VxOpt->output_map(stat_pjc), + 1, 1, stat_at, i_stat_row, + txt_at[i_pjc], i_txt_row[i_pjc]); + } } - } - // Write PRC output - if(conf_info.VxOpt[i_vx].OutputMap[txt_file_type[i_pjc]] != STATOutputType_None) { + // Write PRC output + if(pgi.VxOpt->output_map(stat_pjc) != STATOutputType_None) { - if(conf_info.VxOpt[i_vx].DevFlag) { - write_pct_row(shc, pct_dev, - conf_info.VxOpt[i_vx].OutputMap[txt_file_type[i_prc]], - 1, 1, stat_at, i_stat_row, - txt_at[i_prc], i_txt_row[i_prc]); - } - if(conf_info.VxOpt[i_vx].OpsFlag) { - write_pct_row(shc, pct_ops, - conf_info.VxOpt[i_vx].OutputMap[txt_file_type[i_prc]], - 1, 1, stat_at, i_stat_row, - txt_at[i_prc], i_txt_row[i_prc]); + if(pgi.VxOpt->DevFlag) { + shc.set_fcst_var(probgen_dev_name); + shc.set_obs_var (probgen_dev_name); + write_pct_row(shc, pgi.PCTDev[lead_hr], + pgi.VxOpt->output_map(stat_pjc), + 1, 1, stat_at, i_stat_row, + txt_at[i_prc], i_txt_row[i_prc]); + } + if(pgi.VxOpt->OpsFlag) { + shc.set_fcst_var(probgen_ops_name); + shc.set_obs_var (probgen_ops_name); + write_pct_row(shc, pgi.PCTOps[lead_hr], + pgi.VxOpt->output_map(stat_pjc), + 1, 1, stat_at, i_stat_row, + txt_at[i_prc], i_txt_row[i_prc]); + } } - } + } // end for i return; } diff --git a/met/src/tools/tc_utils/tc_gen/tc_gen.h b/met/src/tools/tc_utils/tc_gen/tc_gen.h index fc2ce357eb..085b2a7cf8 100644 --- a/met/src/tools/tc_utils/tc_gen/tc_gen.h +++ b/met/src/tools/tc_utils/tc_gen/tc_gen.h @@ -80,6 +80,8 @@ static const char *txt_file_abbr[n_txt] = { const ConcatString genesis_name ("GENESIS"); const ConcatString genesis_dev_name("GENESIS_DEV"); const ConcatString genesis_ops_name("GENESIS_OPS"); +const ConcatString probgen_dev_name("PROBGEN_DEV"); +const ConcatString probgen_ops_name("PROBGEN_OPS"); // Maximum Best track cyclone number to be processed // Cyclone numbers > 50 are for testing or invests diff --git a/met/src/tools/tc_utils/tc_gen/tc_gen_conf_info.cc b/met/src/tools/tc_utils/tc_gen/tc_gen_conf_info.cc index ee3b30679a..389bf5d603 100644 --- a/met/src/tools/tc_utils/tc_gen/tc_gen_conf_info.cc +++ b/met/src/tools/tc_utils/tc_gen/tc_gen_conf_info.cc @@ -164,7 +164,7 @@ void TCGenVxOpt::clear() { OpsHitBeg = OpsHitEnd = bad_data_int; DiscardFlag = false; DevFlag = OpsFlag = false; - GenProbThresh.clear(); + ProbGenThresh.clear(); CIAlpha = bad_data_double; OutputMap.clear(); NcInfo.clear(); @@ -290,10 +290,10 @@ void TCGenVxOpt::process_config(Dictionary &dict) { exit(1); } - // Conf: genesis_prob_thresh - GenProbThresh = dict.lookup_thresh_array(conf_key_genesis_prob_thresh); - GenProbThresh = string_to_prob_thresh(GenProbThresh.get_str().c_str()); - check_prob_thresh(GenProbThresh); + // Conf: prob_genesis_thresh + ProbGenThresh = dict.lookup_thresh_array(conf_key_prob_genesis_thresh); + ProbGenThresh = string_to_prob_thresh(ProbGenThresh.get_str().c_str()); + check_prob_thresh(ProbGenThresh); // Conf: ci_alpha CIAlpha = dict.lookup_double(conf_key_ci_alpha); @@ -891,13 +891,12 @@ int TCGenConfInfo::get_max_n_prob_thresh() const { int i, n; for(i=0,n=0; iValidGenesisDHrThresh; + + // Setup the default PCTInfo object + DefaultPCT.set_fthresh(VxOpt->ProbGenThresh); + DefaultPCT.allocate_n_alpha(1); + DefaultPCT.alpha[0] = VxOpt->CIAlpha; + + return; +} + +//////////////////////////////////////////////////////////////////////// + +void ProbGenPCTInfo::add(const ProbGenInfo &pgi) { + int i, lead_hr; + + // Store the model name + if(Model.empty()) Model = pgi.technique(); + + // Track the range of valid times + if(InitBeg == 0 || InitBeg > pgi.init()) InitBeg = pgi.init(); + if(InitEnd == 0 || InitEnd < pgi.init()) InitEnd = pgi.init(); + + // Add new map entry for each lead time, if needed + for(i=0; i OutputMap; TCGenNcOutInfo NcInfo; @@ -315,6 +315,41 @@ class GenCTCInfo { //////////////////////////////////////////////////////////////////////// +class ProbGenPCTInfo { + + private: + + void init_from_scratch(); + + PCTInfo DefaultPCT; + + public: + + ProbGenPCTInfo(); + ~ProbGenPCTInfo(); + + ////////////////////////////////////////////////////////////////// + + ConcatString Model; + unixtime InitBeg, InitEnd; + const TCGenVxOpt* VxOpt; + IntArray LeadTimes; + map PCTDev, PCTOps; + + SingleThresh ValidGenesisDHrThresh; + + ////////////////////////////////////////////////////////////////// + + void clear(); + + void set_vx_opt(const TCGenVxOpt *); + + void add(const ProbGenInfo &); + +}; + +//////////////////////////////////////////////////////////////////////// + #endif /* __TC_GEN_CONF_INFO_H__ */ ////////////////////////////////////////////////////////////////////////