diff --git a/src/tools/other/mode_time_domain/mtd.cc b/src/tools/other/mode_time_domain/mtd.cc index 1a1fcf2f1d..7038fee0d8 100644 --- a/src/tools/other/mode_time_domain/mtd.cc +++ b/src/tools/other/mode_time_domain/mtd.cc @@ -17,9 +17,10 @@ // ---- ---- ---- ----------- // 000 10-15-15 Bullock New // 001 05-15-17 Prestopnik P. Added regrid shape -// 002 04-18-19 Halley Gotway Add FCST and OBS units. -// 003 04-25-19 Halley Gotway Add percentiles to 2D output. -// 004 07-06-22 Howard Soh METplus-Internal #19 Rename main to met_ain +// 002 04-18-19 Halley Gotway Add FCST and OBS units +// 003 04-25-19 Halley Gotway Add percentiles to 2D output +// 004 07-06-22 Howard Soh METplus-Internal #19 Rename main to met_main +// 005 08-01-22 Albo MET #1971 Differing time steps // //////////////////////////////////////////////////////////////////////// @@ -205,13 +206,20 @@ MtdIntFile fcst_obj, obs_obj; MM_Engine engine; + // + // storage for valid times + // + +vector valid_times_fcst; +vector valid_times_obs; + // // read the data files // -mtd_read_data(config, *(config.fcst_info), fcst_filenames, fcst_raw); +valid_times_fcst = mtd_read_data(config, *(config.fcst_info), fcst_filenames, fcst_raw); -mtd_read_data(config, *(config.obs_info), obs_filenames, obs_raw); +valid_times_obs = mtd_read_data(config, *(config.obs_info), obs_filenames, obs_raw); if ( fcst_raw.nt() != obs_raw.nt() ) { @@ -436,13 +444,21 @@ for (j=0; j<(fcst_obj.n_objects()); ++j) { mask_2d = fcst_obj.const_t_mask(t, j + 1); // 1-based + if (t < 0 || t >= (int)valid_times_fcst.size()) { + mlog << Error + << "\n " << program_name + << ": index " << t << " out of forecast valid times range 0 to " + << valid_times_fcst.size()-1 << "\n\n"; + exit ( 1 ); + } + fcst_raw.get_data_plane(t, raw_2d); att_2 = calc_2d_single_atts(mask_2d, raw_2d, j + 1, config.inten_perc_value); att_2.set_fcst(); - att_2.set_valid_time(fcst_obj.start_valid_time() + t*(fcst_obj.delta_t())); + att_2.set_valid_time(valid_times_fcst[t]); att_2.set_lead_time(fcst_obj.lead_time(t)); @@ -470,13 +486,21 @@ for (j=0; j<(obs_obj.n_objects()); ++j) { mask_2d = obs_obj.const_t_mask(t, j + 1); // 1-based + if (t < 0 || t >= (int)valid_times_obs.size()) { + mlog << Error + << "\n " << program_name << ": index " << t + << " out of obs valid times range 0 to " + << valid_times_obs.size()-1 << "\n\n"; + exit ( 1 ); + } + obs_raw.get_data_plane(t, raw_2d); att_2 = calc_2d_single_atts(mask_2d, raw_2d, j + 1, config.inten_perc_value); att_2.set_obs(); - att_2.set_valid_time(obs_obj.start_valid_time() + t*(obs_obj.delta_t())); + att_2.set_valid_time(valid_times_obs[t]); att_2.set_lead_time(obs_obj.lead_time(t)); @@ -680,10 +704,15 @@ if ( have_pairs ) { for (t=(att_3.tmin()); t<=(att_3.tmax()); ++t) { - // mask_2d = mask.const_t_mask(t, j + 1); // 1-based mask_2d = mask.const_t_mask(t, 1); // 1-based - // cout << "j = " << j << ", vol = " << mask_2d.object_volume(0) << '\n'; + if (t < 0 || t >= (int)valid_times_fcst.size()) { + mlog << Error + << "\n " << program_name << ": index " << t + << " out of forecast valid times range 0 to " + << valid_times_fcst.size()-1 << "\n\n"; + exit ( 1 ); + } fcst_raw.get_data_plane(t, raw_2d); @@ -691,7 +720,7 @@ if ( have_pairs ) { att_2.set_fcst(); - att_2.set_valid_time(fcst_obj.start_valid_time() + t*(fcst_obj.delta_t())); + att_2.set_valid_time(valid_times_fcst[t]); att_2.set_lead_time(fcst_obj.lead_time(t)); @@ -723,16 +752,23 @@ if ( have_pairs ) { for (t=(att_3.tmin()); t<=(att_3.tmax()); ++t) { - // mask_2d = mask.const_t_mask(t, j + 1); // 1-based mask_2d = mask.const_t_mask(t, 1); // 1-based + if (t < 0 || t >= (int)valid_times_obs.size()) { + mlog << Error + << "\n " << program_name << ": index " << t + << " out of obs valid times range 0 to " + << valid_times_obs.size()-1 << "\n\n"; + exit ( 1 ); + } + obs_raw.get_data_plane(t, raw_2d); att_2 = calc_2d_single_atts(mask_2d, raw_2d, j + 1, config.inten_perc_value); att_2.set_obs(); - att_2.set_valid_time(obs_obj.start_valid_time() + t*(obs_obj.delta_t())); + att_2.set_valid_time(valid_times_obs[t]); att_2.set_lead_time(obs_obj.lead_time(t)); @@ -785,9 +821,8 @@ mlog << Debug(2) << "Creating 2D constant-time slice attributes file: \"" << path << "\"\n"; -do_2d_txt_output(fcst_raw, obs_raw, - fcst_simple_att_2d, obs_simple_att_2d, - fcst_cluster_att_2d, obs_cluster_att_2d, config, path.c_str()); + do_2d_txt_output(fcst_raw, obs_raw, fcst_simple_att_2d, obs_simple_att_2d, + fcst_cluster_att_2d, obs_cluster_att_2d, config, path.c_str()); // // write simple single attributes @@ -872,7 +907,6 @@ mlog << Debug(2) do_mtd_nc_output(config.nc_info, engine, fcst_raw, obs_raw, fcst_obj, obs_obj, config, path.c_str()); - // // done // @@ -1058,11 +1092,17 @@ ConcatString prefix; ConcatString path; + // + // storage for valid times + // + +vector valid_times; + // // read the data files // -mtd_read_data(config, *(config.fcst_info), single_filenames, raw); +valid_times = mtd_read_data(config, *(config.fcst_info), single_filenames, raw); // // copy forecast name/units/level to observation @@ -1166,13 +1206,21 @@ for (j=0; j<(obj.n_objects()); ++j) { mask_2d = obj.const_t_mask(t, j + 1); // 1-based + if (t < 0 || t >= (int)valid_times.size()) { + mlog << Error + << "\n " << program_name << ": index " << t + << " out of valid times range 0 to " + << valid_times.size()-1 << "\n\n"; + exit ( 1 ); + } + raw.get_data_plane(t, raw_2d); att_2 = calc_2d_single_atts(mask_2d, raw_2d, j + 1, config.inten_perc_value); att_2.set_fcst(); - att_2.set_valid_time(obj.start_valid_time() + t*(obj.delta_t())); + att_2.set_valid_time(valid_times[t]); att_2.set_lead_time(obj.lead_time(t)); @@ -1232,7 +1280,6 @@ mlog << Debug(2) do_mtd_nc_output(config.nc_info, raw, obj, config, path.c_str()); - // // done // diff --git a/src/tools/other/mode_time_domain/mtd_file_base.cc b/src/tools/other/mode_time_domain/mtd_file_base.cc index 4295a2ba85..f1e3571d17 100644 --- a/src/tools/other/mode_time_domain/mtd_file_base.cc +++ b/src/tools/other/mode_time_domain/mtd_file_base.cc @@ -91,6 +91,8 @@ Nx = Ny = Nt = 0; StartValidTime = (unixtime) 0; +ActualValidTimes.clear(); + DeltaT = 0; Filename.clear(); @@ -292,7 +294,6 @@ return; } - //////////////////////////////////////////////////////////////////////// @@ -306,6 +307,19 @@ return; } +//////////////////////////////////////////////////////////////////////// + + +void MtdFileBase::init_actual_valid_times(const vector &validTimes) + +{ + +ActualValidTimes = validTimes; + +return; + +} + //////////////////////////////////////////////////////////////////////// @@ -349,6 +363,25 @@ return ( StartValidTime + t*DeltaT ); } +//////////////////////////////////////////////////////////////////////// + + +unixtime MtdFileBase::actual_valid_time(int t) const + +{ + +if ( (t < 0) || ( t >= (int)ActualValidTimes.size()) ) { + + mlog << Error << "\n\n MtdFileBase::valid_time(int t) -> range check error\n\n"; + + exit ( 1 ); + +} + +return ( ActualValidTimes[t] ); + +} + //////////////////////////////////////////////////////////////////////// diff --git a/src/tools/other/mode_time_domain/mtd_file_base.h b/src/tools/other/mode_time_domain/mtd_file_base.h index 205b5feb92..98d30f10a1 100644 --- a/src/tools/other/mode_time_domain/mtd_file_base.h +++ b/src/tools/other/mode_time_domain/mtd_file_base.h @@ -85,9 +85,11 @@ class MtdFileBase { int Nx, Ny, Nt; - unixtime StartValidTime; + unixtime StartValidTime; // useful for constant time increments - int DeltaT; // seconds + int DeltaT; // seconds, useful for constant time increments + + vector ActualValidTimes; // useful for uneven time increments IntArray Lead_Times; @@ -113,8 +115,11 @@ class MtdFileBase { void set_grid(const Grid &); void set_start_valid_time (unixtime); + void set_delta_t (int); // seconds + void init_actual_valid_times(const vector &validTimes); + void set_lead_time(int index, int value); void set_filetype(MtdFileType); @@ -143,6 +148,8 @@ class MtdFileBase { unixtime valid_time (int) const; + unixtime actual_valid_time (int) const; + int lead_time (int index) const; // diff --git a/src/tools/other/mode_time_domain/mtd_read_data.cc b/src/tools/other/mode_time_domain/mtd_read_data.cc index 2b4eb852d0..d40335594b 100644 --- a/src/tools/other/mode_time_domain/mtd_read_data.cc +++ b/src/tools/other/mode_time_domain/mtd_read_data.cc @@ -22,12 +22,13 @@ using namespace std; #include "mtd_read_data.h" +#include "num_array.h" //////////////////////////////////////////////////////////////////////// -void mtd_read_data(MtdConfigInfo & config, VarInfo & varinfo, - const StringArray & filenames, MtdFloatFile & raw) +vector mtd_read_data(MtdConfigInfo & config, VarInfo & varinfo, + const StringArray & filenames, MtdFloatFile & raw) { @@ -39,16 +40,11 @@ if ( filenames.n() < 2 ) { } -int j; +int j, k; Met2dDataFile * data_2d_file = 0; Met2dDataFileFactory factory; DataPlane plane; -unixtime * valid_times = 0; - - - -valid_times = new unixtime [filenames.n()]; - +vector valid_times; // // read the files // @@ -76,7 +72,7 @@ for (j=0; j<(filenames.n()); ++j) { } - valid_times[j] = plane.valid(); + valid_times.push_back(plane.valid()); if ( j == 0 ) { @@ -96,30 +92,76 @@ for (j=0; j<(filenames.n()); ++j) { } // for j -varinfo.set_valid(valid_times[0]); + varinfo.set_valid(valid_times[0]); + + // store the valid times vector into the raw data for later use in do_2d_txt_output() + + raw.init_actual_valid_times(valid_times); // - // check the time intervals + // check the time intervals for consistency + // Store the time differences between succesive valid times in an array + // See if differences are constant or not, and if not see if all diffs are months // -unixtime dt_start, dt; +unixtime dt_start; +vector dtVector; dt_start = valid_times[1] - valid_times[0]; +dtVector.push_back(dt_start); -for (j=2; j<(filenames.n()); ++j) { - - dt = valid_times[j] - valid_times[j - 1]; - - if ( dt != dt_start ) { - - mlog << Error << "\n\n mtd_read_data() -> file time increments are not constant!\n\n"; - - exit ( 1 ); +for (size_t k=2; k File time increments are not constant, could be problematic\n"; + mlog << Warning << " mtd_read_data() -> Using MODE of the increments, mode=" << dt_start << "\n"; + mlog << Warning << " mtd_read_data() -> Time increment properties: mean=" << umean << " variance=" << uvar << " sqrt(var)=" << suvar << "\n\n"; + } } - + // // load up the rest of the MtdFloatFile class members // @@ -138,9 +180,7 @@ raw.calc_data_minmax(); // done // -if ( valid_times ) { delete [] valid_times; valid_times = 0; } - -return; +return valid_times; } diff --git a/src/tools/other/mode_time_domain/mtd_read_data.h b/src/tools/other/mode_time_domain/mtd_read_data.h index 26d0a9ba18..0da67cc2a9 100644 --- a/src/tools/other/mode_time_domain/mtd_read_data.h +++ b/src/tools/other/mode_time_domain/mtd_read_data.h @@ -26,9 +26,13 @@ //////////////////////////////////////////////////////////////////////// +// +// returns the actual valid times +// -extern void mtd_read_data(MtdConfigInfo &, VarInfo &, - const StringArray & filenames, MtdFloatFile &); +extern vector mtd_read_data(MtdConfigInfo &, VarInfo &, + const StringArray & filenames, + MtdFloatFile &); //////////////////////////////////////////////////////////////////////// diff --git a/src/tools/other/mode_time_domain/mtd_txt_output.cc b/src/tools/other/mode_time_domain/mtd_txt_output.cc index a0bf061117..2372823d25 100644 --- a/src/tools/other/mode_time_domain/mtd_txt_output.cc +++ b/src/tools/other/mode_time_domain/mtd_txt_output.cc @@ -496,7 +496,7 @@ for (j=0; j<(fcst_simple_att.n()); ++j) { table.set_entry(r, fcst_valid_column, unix_to_yyyymmdd_hhmmss(fcst_simple_att.valid_time(j))); - table.set_entry(r, obs_valid_column, unix_to_yyyymmdd_hhmmss(obs_raw.valid_time(t))); + table.set_entry(r, obs_valid_column, unix_to_yyyymmdd_hhmmss(obs_raw.actual_valid_time(t))); table.set_entry(r, obs_lead_column, sec_to_hhmmss(obs_raw.lead_time(t))); @@ -508,7 +508,7 @@ for (j=0; j<(obs_simple_att.n()); ++j) { t = obs_simple_att.time_index(j); - table.set_entry(r, fcst_valid_column, unix_to_yyyymmdd_hhmmss(fcst_raw.valid_time(t))); + table.set_entry(r, fcst_valid_column, unix_to_yyyymmdd_hhmmss(fcst_raw.actual_valid_time(t))); table.set_entry(r, fcst_lead_column, sec_to_hhmmss(fcst_raw.lead_time(t))); @@ -528,7 +528,7 @@ for (j=0; j<(fcst_cluster_att.n()); ++j) { table.set_entry(r, fcst_valid_column, unix_to_yyyymmdd_hhmmss(fcst_cluster_att.valid_time(j))); - table.set_entry(r, obs_valid_column, unix_to_yyyymmdd_hhmmss(obs_raw.valid_time(t))); + table.set_entry(r, obs_valid_column, unix_to_yyyymmdd_hhmmss(obs_raw.actual_valid_time(t))); table.set_entry(r, obs_lead_column, sec_to_hhmmss(obs_raw.lead_time(t))); @@ -540,7 +540,7 @@ for (j=0; j<(obs_cluster_att.n()); ++j) { t = obs_cluster_att.time_index(j); - table.set_entry(r, fcst_valid_column, unix_to_yyyymmdd_hhmmss(fcst_raw.valid_time(t))); + table.set_entry(r, fcst_valid_column, unix_to_yyyymmdd_hhmmss(fcst_raw.actual_valid_time(t))); table.set_entry(r, fcst_lead_column, sec_to_hhmmss(fcst_raw.lead_time(t)));