Bugfix #2412 main_v11.0 climo

internal/scripts/docker/build_met_docker.sh

@@ -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

src/basic/vx_cal/doyhms_to_unix.cc

@@ -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;
@@ -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 );
+
+}
+
+
 ////////////////////////////////////////////////////////////////////////

src/basic/vx_cal/vx_cal.h

@@ -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
@@ -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);

src/basic/vx_util/interp_util.cc

@@ -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 "

src/libcode/vx_statistics/read_climo.cc

@@ -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))) {
@@ -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
 
@@ -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]),