Skip to content

Commit

Permalink
Bugfix #2412 main_v11.0 climo (#2420)
Browse files Browse the repository at this point in the history 
Co-authored-by: MET Tools Test Account <met_test@seneca.rap.ucar.edu>
  • Loading branch information
@JohnHalleyGotway
JohnHalleyGotway and MET Tools Test Account authored Jan 25, 2023
1 parent c3b3b9a commit 5d36e95
Show file tree
Hide file tree
Showing 5 changed files with 127 additions and 81 deletions.
4 changes: 2 additions & 2 deletions internal/scripts/docker/build_met_docker.sh
View file
Original file line number Diff line number Diff line change
Expand Up @@ -15,14 +15,14 @@ fi

LOG_FILE=/met/MET-${MET_GIT_NAME}/make_install.log
echo "Compiling MET ${MET_GIT_NAME} and writing log file ${LOG_FILE}"
make install > ${LOG_FILE}
make -j install > ${LOG_FILE}
if [ $? != 0 ]; then
exit 1
fi

LOG_FILE=/met/logs/MET-${MET_GIT_NAME}_make_test.log
echo "Testing MET ${MET_GIT_NAME} and writing log file ${LOG_FILE}"
make test > ${LOG_FILE} 2>&1
make -j test > ${LOG_FILE} 2>&1
if [ $? != 0 ]; then
exit 1
fi
Expand Down
65 changes: 65 additions & 0 deletions src/basic/vx_cal/doyhms_to_unix.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -97,8 +97,28 @@ return ( s );

}


////////////////////////////////////////////////////////////////////////


int unix_to_day_of_year(unixtime u) {

int mon, day, yr, hr, min, sec;

unix_to_mdyhms(u, mon, day, yr, hr, min, sec);

int sec_of_year = mdyhms_to_unix(mon, day, 1970, 0, 0, 0);

int sec_per_day = 60 * 60 * 24;

return ( sec_of_year / sec_per_day );

}


////////////////////////////////////////////////////////////////////////


long unix_to_long_yyyymmddhh(unixtime u) {

int mon, day, yr, hr, min, sec;
Expand All @@ -112,4 +132,49 @@ return ( yyyymmddhh );

}


////////////////////////////////////////////////////////////////////////


int sec_of_day_diff(unixtime ut1, unixtime ut2) {

int sec_per_day = 60 * 60 * 24;

int s1 = unix_to_sec_of_day(ut1);

int s2 = unix_to_sec_of_day(ut2);

int dt = s2 - s1;

if ( dt < -1 * sec_per_day/2 ) dt += sec_per_day;
else if ( dt > sec_per_day/2 ) dt -= sec_per_day;

return ( dt );

}


////////////////////////////////////////////////////////////////////////


int day_of_year_diff(unixtime ut1, unixtime ut2) {


int sec_per_day = 60 * 60 * 24;
int day_per_year = 365;

int d1 = unix_to_day_of_year(ut1);

int d2 = unix_to_day_of_year(ut2);

int dt = d2 - d1;

if ( dt < -1 * day_per_year/2.0 ) dt += day_per_year;
else if ( dt > day_per_year/2.0 ) dt -= day_per_year;

return ( dt );

}


////////////////////////////////////////////////////////////////////////
10 changes: 10 additions & 0 deletions src/basic/vx_cal/vx_cal.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -81,6 +81,8 @@ extern int hms_to_sec (int hour, int min, int sec);

extern int unix_to_sec_of_day (unixtime u);

extern int unix_to_day_of_year (unixtime u);

extern long unix_to_long_yyyymmddhh (unixtime u);

// Parse time strings
Expand Down Expand Up @@ -128,6 +130,14 @@ extern ConcatString sec_to_timestring(int);
////////////////////////////////////////////////////////////////////////


extern int sec_of_day_diff (unixtime ut1, unixtime ut2);

extern int day_of_year_diff (unixtime ut1, unixtime ut2);


////////////////////////////////////////////////////////////////////////


extern bool is_datestring(const char * text);

extern bool is_yyyymmdd(const char * text);
Expand Down
10 changes: 6 additions & 4 deletions src/basic/vx_util/interp_util.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1280,14 +1280,16 @@ DataPlane valid_time_interp(const DataPlane &in1, const DataPlane &in2,
double w1 = bad_data_double;
double w2 = bad_data_double;

// Store min and max valid times.
// Store min and max valid times
dp1 = (in1.valid() <= in2.valid() ? in1 : in2);
dp2 = (in1.valid() > in2.valid() ? in1 : in2);

// Check for matching valid times
if(dp1.valid() == dp2.valid()) return(dp1);

// Range check the times
if(dp1.valid() > to_ut || dp2.valid() < to_ut ||
dp1.valid() == dp2.valid()) {
mlog << Error << "\time_interp() -> "
if(dp1.valid() > to_ut || dp2.valid() < to_ut) {
mlog << Error << "\ntime_interp() -> "
<< "the interpolation time " << unix_to_yyyymmdd_hhmmss(to_ut)
<< " must fall between the input times: "
<< unix_to_yyyymmdd_hhmmss(dp1.valid()) << " and "
Expand Down
119 changes: 44 additions & 75 deletions src/libcode/vx_statistics/read_climo.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -150,21 +150,19 @@ void read_climo_file(const char *climo_file, GrdFileType ctype,
int day_ts, int hour_ts, const Grid &vx_grid,
const RegridInfo &regrid_info,
DataPlaneArray &dpa) {

Met2dDataFileFactory mtddf_factory;
Met2dDataFile *mtddf = (Met2dDataFile *) 0;

VarInfoFactory info_factory;
VarInfo *info = (VarInfo *) 0;

DataPlaneArray cur_dpa;
DataPlaneArray clm_dpa;
DataPlane dp;

int n_climo, i;
int vld_mon, vld_day, vld_yr, vld_hr, vld_min, vld_sec;
int clm_mon, clm_day, clm_yr, clm_hr, clm_min, clm_sec;
unixtime ut;
int i, n_clm, day_diff_sec, hms_diff_sec;

ConcatString cur_ut_cs;
ConcatString clm_ut_cs;

// Allocate memory for data file
if(!(mtddf = mtddf_factory.new_met_2d_data_file(climo_file, ctype))) {
Expand All @@ -179,94 +177,65 @@ void read_climo_file(const char *climo_file, GrdFileType ctype,
info->set_dict(*dict);

// Read data planes
n_climo = mtddf->data_plane_array(*info, cur_dpa);

// Compute the valid month and day
unix_to_mdyhms(vld_ut, vld_mon, vld_day, vld_yr,
vld_hr, vld_min, vld_sec);
n_clm = mtddf->data_plane_array(*info, clm_dpa);

// Loop through matching records
for(i=0; i<n_climo; i++) {

// Store current valid time string
cur_ut_cs = unix_to_yyyymmdd_hhmmss(cur_dpa[i].valid());

// Compute the current month and day
unix_to_mdyhms(cur_dpa[i].valid(), clm_mon, clm_day, clm_yr,
clm_hr, clm_min, clm_sec);

// Recompute the unixtime to check the hour of the day.
// Use the valid YYYYMMDD and climo HHMMSS.
ut = mdyhms_to_unix(vld_mon, vld_day, vld_yr,
clm_hr, clm_min, clm_sec);

// Check the hour time step.
if(!is_bad_data(hour_ts) && labs(ut - vld_ut) >= hour_ts) {
mlog << Debug(3) << "Skipping the " << cur_ut_cs << " \""
<< info->magic_str()
<< "\" climatology field since the time offset ("
<< labs(ut - vld_ut)
<< " seconds) >= the \"" << conf_key_hour_interval
<< "\" entry (" << hour_ts << " seconds) from file: "
<< climo_file << "\n";
for(i=0; i<n_clm; i++) {

// Store climo time string
clm_ut_cs = unix_to_yyyymmdd_hhmmss(clm_dpa[i].valid());

// Compute day and hour time offsets in seconds
day_diff_sec = day_of_year_diff(vld_ut, clm_dpa[i].valid()) * sec_per_day;
hms_diff_sec = sec_of_day_diff (vld_ut, clm_dpa[i].valid());

// Check the day time step
if(!is_bad_data(day_ts) && abs(day_diff_sec) >= day_ts) {
mlog << Debug(3) << "Skipping " << clm_ut_cs << " \"" << info->magic_str()
<< "\" climatology field with " << day_diff_sec / sec_per_day
<< " day offset (" << conf_key_day_interval << " = "
<< day_ts / sec_per_day << ") from file \""
<< climo_file << "\".\n";
continue;
}

// Recompute the unixtime to check the day of the year.
// Use the valid YYYY and climo MMDD_HHMMSS.
ut = mdyhms_to_unix(clm_mon, clm_day, vld_yr,
clm_hr, clm_min, clm_sec);

// Check the day time step.
if(!is_bad_data(day_ts)) {

// For daily climatology, check the hour timestep.
if(day_ts <= 3600*24 && labs(ut - vld_ut) >= hour_ts) {
mlog << Debug(3) << "Skipping the " << cur_ut_cs << " \""
<< info->magic_str()
<< "\" climatology field since the time offset ("
<< labs(ut - vld_ut) << " seconds) >= the \""
<< conf_key_hour_interval << "\" entry (" << hour_ts
<< " seconds) from daily climatology file: "
<< climo_file << "\n";
continue;
}
// For non-daily climatology, check the day timestep.
else if(labs(ut - vld_ut) >= day_ts) {
mlog << Debug(3) << "Skipping the " << cur_ut_cs << " \""
<< info->magic_str()
<< "\" climatology field since the time offset ("
<< labs(ut - vld_ut) << " seconds) >= the \""
<< conf_key_day_interval << "\" entry (" << day_ts
<< " seconds) from file: " << climo_file << "\n";
continue;
}
// Check the hour time step
if(!is_bad_data(hour_ts) && abs(hms_diff_sec) >= hour_ts) {
mlog << Debug(3) << "Skipping " << clm_ut_cs << " \"" << info->magic_str()
<< "\" climatology field with " << (double) hms_diff_sec / sec_per_hour
<< " hour offset (" << conf_key_hour_interval << " = "
<< hour_ts / sec_per_hour << ") from file \""
<< climo_file << "\".\n";
continue;
}

// Compute climo timestamp relative to the valid time
unixtime clm_vld_ut = vld_ut + day_diff_sec + hms_diff_sec;

// Print log message for matching record
mlog << Debug(4)
<< "Found matching " << cur_ut_cs << " \""
<< info->magic_str() << "\" climatology field in file \""
<< climo_file << "\".\n";
mlog << Debug(3) << "Storing " << clm_ut_cs << " \"" << info->magic_str()
<< "\" climatology field with " << day_diff_sec / sec_per_day
<< " day, " << (double) hms_diff_sec / sec_per_hour << " hour offset as time "
<< unix_to_yyyymmdd_hhmmss(clm_vld_ut) << " from file \""
<< climo_file << "\".\n";

// Regrid, if needed
if(!(mtddf->grid() == vx_grid)) {
mlog << Debug(2)
<< "Regridding the " << cur_ut_cs << " \""
mlog << Debug(2) << "Regridding " << clm_ut_cs << " \""
<< info->magic_str()
<< "\" climatology field to the verification grid.\n";
dp = met_regrid(cur_dpa[i], mtddf->grid(), vx_grid,
dp = met_regrid(clm_dpa[i], mtddf->grid(), vx_grid,
regrid_info);
}
else {
dp = cur_dpa[i];
dp = clm_dpa[i];
}

// Set the climo time as valid YYYY and climo MMDD_HHMMSS.
dp.set_valid(ut);
// Set the climo time relative to the valid time
dp.set_valid(clm_vld_ut);

// Store the match
dpa.add(dp, cur_dpa.lower(i), cur_dpa.upper(i));
dpa.add(dp, clm_dpa.lower(i), clm_dpa.upper(i));

} // end for i

Expand Down Expand Up @@ -376,7 +345,7 @@ DataPlaneArray climo_time_interp(const DataPlaneArray &dpa, int day_ts,

// For equality, do a single time interpolation.
if(prv_hms == nxt_hms) {
ut = (vld_ut / sec_per_day) + prv_hms;
ut = (vld_ut / sec_per_day)*sec_per_day + prv_hms;
interp_dpa.add(climo_hms_interp(
dpa, it->second, ut, mthd),
dpa.lower(it->second[0]),
Expand Down

0 comments on commit 5d36e95

Please sign in to comment.