Merge pull request #212 from UW-Hydro/develop

Prep for 2.2.1 release

docs/whats-new.rst

@@ -5,11 +5,20 @@ What's New
 
 .. _whats-new.2.2.0:
 
+v2.2.l
+------
+Bug fixes
+~~~~~~~~~
+- Fixed bug where timestamps got duplicated for sub-hourly
+  disaggregation time periods
+- Fixed bug where polar latitudes caused error in setting
+  the rise and set times for temperature disaggregation
+
 v2.2.0
 -------
 Enhancements
 ~~~~~~~~~~~~
-- Can now specify ``period_ending`` in the configuration to move 
+- Can now specify ``period_ending`` in the configuration to move
   timestamps to end of period instead of beginning of period
 - Addition of tutorial in README.md and main documentation
 - Addition of paper to be submitted to JOSS

metsim/disaggregate.py

@@ -28,6 +28,7 @@
 from metsim.datetime import date_range
 from metsim.physics import svp
 
+import sys
 
 def disaggregate(df_daily: pd.DataFrame, params: dict,
                  solar_geom: dict, t_begin: list=None,
@@ -77,6 +78,7 @@ def disaggregate(df_daily: pd.DataFrame, params: dict,
 
     t_Tmin, t_Tmax = set_min_max_hour(solar_geom['tiny_rad_fract'],
                                       df_daily.index.dayofyear - 1,
+                                      solar_geom['daylength'],
                                       n_days, ts, params)
 
     df_disagg['temp'] = temp(
@@ -114,8 +116,8 @@ def disaggregate(df_daily: pd.DataFrame, params: dict,
     return df_disagg.fillna(method='ffill').fillna(method='bfill')
 
 
-def set_min_max_hour(tiny_rad_fract: np.array, yday: np.array, n_days: int,
-                     ts: int, params: dict) -> Tuple[np.array]:
+def set_min_max_hour(tiny_rad_fract: np.ndarray, yday: np.ndarray, daylength: np.ndarray,
+                     n_days: int, ts: int, params: dict) -> Tuple[np.ndarray]:
     """
     Determine the time at which min and max temp
     is reached for each day.
@@ -146,11 +148,26 @@ def set_min_max_hour(tiny_rad_fract: np.array, yday: np.array, n_days: int,
     # calculate minute of sunrise and sunset for each day of the year
     rad_mask = 1 * (tiny_rad_fract > 0)
     mask = np.diff(rad_mask)
-    rise_times = np.where(mask > 0)[1] * (cnst.SW_RAD_DT / cnst.SEC_PER_MIN)
-    set_times = np.where(mask < 0)[1] * (cnst.SW_RAD_DT / cnst.SEC_PER_MIN)
+
+    # north of the polar circle radiation values of mask </> 0 are eleminated for
+    # sunset/sunrise resulting in an array containing less than 365 days for one year
+    rise_times = np.zeros(mask.shape[0])
+    loc, mult = np.where(mask > 0)
+    for i, j in zip(loc, mult):
+        rise_times[i] = j * (cnst.SW_RAD_DT / cnst.SEC_PER_MIN)
+
+    set_times = np.zeros(mask.shape[0])
+    loc, mult = np.where(mask < 0)
+    for i, j in zip(loc, mult):
+        set_times[i] = j * (cnst.SW_RAD_DT / cnst.SEC_PER_MIN)
+
+    # Set rise and set times to be equally spaced when in polar region
+    day_tot = ((cnst.SEC_PER_DAY / cnst.SW_RAD_DT)
+               * (cnst.SW_RAD_DT / cnst.SEC_PER_MIN))
+    rise_times[rise_times == 0] = day_tot / 6
+    set_times[set_times == 0] = 4 * day_tot / 6
 
     if params['utc_offset']:
-        # not used elsewhere:
         # rad_fract_per_day = int(cnst.SEC_PER_DAY/cnst.SW_RAD_DT)
         utc_offset = int(((params.get("lon", 0) - params.get("theta_s", 0)) /
                           cnst.DEG_PER_REV) * cnst.MIN_PER_DAY)

metsim/metsim.py

@@ -70,31 +70,44 @@
 
 attrs = {'pet': {'units': 'mm timestep-1',
                  'long_name': 'potential evaporation',
-                 'standard_name': 'water_potential_evaporation_flux'},
+                 'standard_name': 'water_potential_evaporation_flux',
+                 'missing_value': np.nan, 'fill_value': np.nan},
          'prec': {'units': 'mm timestep-1', 'long_name': 'precipitation',
-                  'standard_name': 'precipitation_flux'},
+                  'standard_name': 'precipitation_flux',
+                  'missing_value': np.nan, 'fill_value': np.nan},
          'shortwave': {'units': 'W m-2', 'long_name': 'shortwave radiation',
-                       'standard_name': 'surface_downwelling_shortwave_flux'},
+                       'standard_name': 'surface_downwelling_shortwave_flux',
+                       'missing_value': np.nan, 'fill_value': np.nan},
          'longwave': {'units': 'W m-2', 'long_name': 'longwave radiation',
-                      'standard_name': 'surface_downwelling_longwave_flux'},
+                      'standard_name': 'surface_downwelling_longwave_flux',
+                      'missing_value': np.nan, 'fill_value': np.nan},
          't_max': {'units': 'C', 'long_name': 'maximum daily air temperature',
-                   'standard_name': 'daily_maximum_air_temperature'},
+                   'standard_name': 'daily_maximum_air_temperature',
+                   'missing_value': np.nan, 'fill_value': np.nan},
          't_min': {'units': 'C', 'long_name': 'minimum daily air temperature',
-                   'standard_name': 'daily_minimum_air_temperature'},
+                   'standard_name': 'daily_minimum_air_temperature',
+                   'missing_value': np.nan, 'fill_value': np.nan},
          'temp': {'units': 'C', 'long_name': 'air temperature',
-                  'standard_name': 'air_temperature'},
+                  'standard_name': 'air_temperature',
+                  'missing_value': np.nan, 'fill_value': np.nan},
          'vapor_pressure': {'units': 'kPa', 'long_name': 'vapor pressure',
-                            'standard_name': 'vapor_pressure'},
+                            'standard_name': 'vapor_pressure',
+                            'missing_value': np.nan, 'fill_value': np.nan},
          'air_pressure': {'units': 'kPa', 'long_name': 'air pressure',
-                          'standard_name': 'air_pressure'},
+                          'standard_name': 'air_pressure',
+                          'missing_value': np.nan, 'fill_value': np.nan},
          'tskc': {'units': 'fraction', 'long_name': 'cloud fraction',
-                  'standard_name': 'cloud_fraction'},
+                  'standard_name': 'cloud_fraction',
+                  'missing_value': np.nan, 'fill_value': np.nan},
          'rel_humid': {'units': '%', 'long_name': 'relative humidity',
-                       'standard_name': 'relative_humidity'},
+                       'standard_name': 'relative_humidity',
+                       'missing_value': np.nan, 'fill_value': np.nan},
          'spec_humid': {'units': '', 'long_name': 'specific humidity',
-                        'standard_name': 'specific_humidity'},
+                        'standard_name': 'specific_humidity',
+                        'missing_value': np.nan, 'fill_value': np.nan},
          'wind': {'units': 'm/s', 'long_name': 'wind speed',
-                  'standard_name': 'wind_speed'},
+                  'standard_name': 'wind_speed',
+                  'missing_value': np.nan, 'fill_value': np.nan},
          '_global': {'conventions': '1.6', 'title': 'Output from MetSim',
                      'institution': 'University of Washington',
                      'source': 'metsim.py',
@@ -408,6 +421,14 @@ def run_slice(self):
         self.setup_output()
         times = self.met_data['time']
         params = self.params.copy()
+        # transform input parameters to floating point values
+        params['sw_prec_thresh'] = float(params['sw_prec_thresh'])
+        params['rain_scalar'] = float(params['rain_scalar'])
+        params['tdew_tol'] = float(params['tdew_tol'])
+        params['tmax_daylength_fraction'] = float(params['tmax_daylength_fraction'])
+        params['tday_coef'] = float(params['tday_coef'])
+        params['tmax_daylength_fraction'] = float(params['tmax_daylength_fraction'])
+        params['lapse_rate'] = float(params['lapse_rate'])
         for index, mask_val in np.ndenumerate(self.domain['mask'].values):
             if mask_val > 0:
                 locs = {d: i for d, i in zip(self.domain['mask'].dims, index)}