From 5ed463d625c3e73df2fb7955be8bce6e2c18eed1 Mon Sep 17 00:00:00 2001 From: grantbuster Date: Fri, 3 Mar 2023 13:21:58 -0700 Subject: [PATCH 1/4] added bias correction method to sam resource handler with tests --- rex/sam_resource.py | 76 +++++++++++++++++++++++++++++++++++++ tests/test_sam_resource.py | 78 ++++++++++++++++++++++++++++++++++++++ 2 files changed, 154 insertions(+) diff --git a/rex/sam_resource.py b/rex/sam_resource.py index a91d8f2b..628d9210 100644 --- a/rex/sam_resource.py +++ b/rex/sam_resource.py @@ -701,6 +701,82 @@ def append_var_list(self, var): self.var_list.append(var) + def bias_correct(self, bc_df): + """Bias correct wind or irradiance data using a table of linear + correction factors per resource gid. + + Parameters + ---------- + bc_df : None | pd.DataFrame + None if not provided or extracted DataFrame with wind or solar + resource bias correction table. This has columns: gid, adder, + scalar. If either adder or scalar is not present, scalar defaults + to 1 and adder to 0. Only windspeed, GHI, and DNI are corrected. + GHI and DNI are corrected with the same correction factors. + """ + + if 'adder' not in bc_df: + logger.info('Bias correction table provided, but "adder" not ' + 'found, defaulting to 0.') + bc_df['adder'] = 0 + + if 'scalar' not in bc_df: + logger.info('Bias correction table provided, but "scalar" not ' + 'found, defaulting to 1.') + bc_df['scalar'] = 1 + + if bc_df.index.name != 'gid': + msg = ('Bias correction table must have "gid" column but only ' + 'found: {}'.format(list(bc_df.columns))) + assert 'gid' in bc_df, msg + bc_df = bc_df.set_index('gid') + + site_arr = np.array(self.sites) + isin = np.isin(self.sites, bc_df.index.values) + if not isin.all(): + missing = site_arr[isin] + msg = ('{} sites were missing from the bias correction table, ' + 'not bias correcting: {}'.format(len(missing), missing)) + logger.warning(msg) + warn(msg) + + scalar = np.expand_dims(bc_df.loc[site_arr[isin], 'scalar'].values, 0) + adder = np.expand_dims(bc_df.loc[site_arr[isin], 'adder'].values, 0) + + if 'windspeed' in self._res_arrays: + logger.debug('Bias correcting windspeeds.') + arr = self._res_arrays['windspeed'] + arr[:, isin] = arr[:, isin] * scalar + adder + arr = np.maximum(arr, 0) + self._res_arrays['windspeed'] = arr + + elif 'ghi' in self._res_arrays: + ghi = self._res_arrays['ghi'] + dni = self._res_arrays['dni'] + dhi = self._res_arrays['dhi'] + ghi_zeros = ghi == 0 + dni_zeros = dni == 0 + dhi_zeros = dhi == 0 + with np.errstate(divide='ignore', invalid='ignore'): + cos_sza = (ghi - dhi) / dni + + ghi[:, isin] = ghi[:, isin] * scalar + adder + dni[:, isin] = dni[:, isin] * scalar + adder + ghi = np.maximum(0, ghi) + dni = np.maximum(0, dni) + ghi[ghi_zeros] = 0 + dni[dni_zeros] = 0 + with np.errstate(divide='ignore', invalid='ignore'): + dhi[:, isin] = ghi[:, isin] - (dni[:, isin] * cos_sza[:, isin]) + + dhi[dni_zeros] = ghi[dni_zeros] + dhi = np.maximum(0, dhi) + dhi[dhi_zeros] = 0 + assert not np.isnan(dhi).any() + self._res_arrays['ghi'] = ghi + self._res_arrays['dni'] = dni + self._res_arrays['dhi'] = dhi + def _check_physical_ranges(self, var, arr, var_slice): """Check physical range of array and enforce usable SAM data. diff --git a/tests/test_sam_resource.py b/tests/test_sam_resource.py index d6ac1100..fa757ede 100644 --- a/tests/test_sam_resource.py +++ b/tests/test_sam_resource.py @@ -410,5 +410,83 @@ def execute_pytest(capture='all', flags='-rapP'): pytest.main(['-q', '--show-capture={}'.format(capture), fname, flags]) +def test_bias_correct_wind(): + """Test linear bias correction functionality on windspeed""" + h5 = os.path.join(TESTDATADIR, 'wtk/ri_100_wtk_2012.h5') + sites = slice(0, 20) + hub_heights = 80 + base_res = WindResource.preload_SAM(h5, sites, hub_heights) + + n = 10 + bc = pd.DataFrame({'gid': np.arange(n), + 'adder': np.random.uniform(-1, 1, n), + 'scalar': np.random.uniform(0.9, 1.1, n)}) + + with pytest.warns() as record: + res = WindResource.preload_SAM(h5, sites, hub_heights) + res.bias_correct(bc) + + assert len(record) == 1 + assert 'missing from the bias correction' in str(record[0].message) + assert np.allclose(res._res_arrays['windspeed'][:, 10:], + base_res._res_arrays['windspeed'][:, 10:]) + assert not (res._res_arrays['windspeed'][:, :10] == + base_res._res_arrays['windspeed'][:, :10]).any() + assert (res._res_arrays['windspeed'] >= 0).all() + + n = 200 + bc = pd.DataFrame({'gid': np.arange(n), + 'adder': np.random.uniform(-1, 1, n), + 'scalar': np.random.uniform(0.9, 1.1, n)}) + + with pytest.warns(None) as record: + res = WindResource.preload_SAM(h5, sites, hub_heights) + res.bias_correct(bc) + + assert not any(record) + assert not (res._res_arrays['windspeed'] == + base_res._res_arrays['windspeed']).any() + assert (res._res_arrays['windspeed'] >= 0).all() + + +def test_bias_correct_solar(): + """Test adder bias correction functionality on irradiance""" + h5 = os.path.join(TESTDATADIR, 'nsrdb/ri_100_nsrdb_2012.h5') + sites = slice(0, 10) + base_res = NSRDB.preload_SAM(h5, sites) + + n = 10 + bc = pd.DataFrame({'gid': np.arange(n), + 'adder': np.random.uniform(-100, 100, n), + 'scalar': np.random.uniform(1, 1, n)}) + + res = NSRDB.preload_SAM(h5, sites) + res.bias_correct(bc) + + for gid in res.sites: + adder = bc.at[gid, 'adder'] + base_ghi = base_res._res_arrays['ghi'][:, gid] + base_dni = base_res._res_arrays['dni'][:, gid] + base_dhi = base_res._res_arrays['dhi'][:, gid] + ghi = res._res_arrays['ghi'][:, gid] + dni = res._res_arrays['dni'][:, gid] + dhi = res._res_arrays['dhi'][:, gid] + assert (ghi >= 0).all() + assert (dni >= 0).all() + assert (dhi >= 0).all() + ghi_mask = (ghi > np.abs(adder)) & (base_ghi > np.abs(adder)) + dni_mask = (dni > np.abs(adder)) & (base_dni > np.abs(adder)) + assert np.allclose(ghi[ghi_mask], base_ghi[ghi_mask] + adder) + assert np.allclose(dni[dni_mask], base_dni[dni_mask] + adder) + + ghi_mask = (ghi > np.abs(adder)) & (base_ghi > np.abs(adder)) + dni_mask = (dni > np.abs(adder)) & (base_dni > np.abs(adder)) + dhi_mask = (dhi > np.abs(adder)) & (base_dhi > np.abs(adder)) + mask = ghi_mask & dni_mask & dhi_mask + cos_sza = (ghi[mask] - dhi[mask]) / (dni[mask]) + base_cos_sza = (base_ghi[mask] - base_dhi[mask]) / (base_dni[mask]) + assert np.allclose(cos_sza, base_cos_sza, atol=0.005) + + if __name__ == '__main__': execute_pytest() From 0c5e76d9f67abe9979dbec571c7d8e843c64accf Mon Sep 17 00:00:00 2001 From: grantbuster Date: Fri, 3 Mar 2023 13:48:30 -0700 Subject: [PATCH 2/4] added recalc of mean resource with bias correctin --- rex/sam_resource.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/rex/sam_resource.py b/rex/sam_resource.py index 628d9210..f2e2c6d8 100644 --- a/rex/sam_resource.py +++ b/rex/sam_resource.py @@ -777,6 +777,10 @@ def bias_correct(self, bc_df): self._res_arrays['dni'] = dni self._res_arrays['dhi'] = dhi + if self._mean_arrays is not None: + for var in self._mean_arrays.keys(): + self._mean_arrays[var] = self._res_arrays[var].mean(axis=0) + def _check_physical_ranges(self, var, arr, var_slice): """Check physical range of array and enforce usable SAM data. From 92481505721aad5a38236e8d8b179fd543c88a8e Mon Sep 17 00:00:00 2001 From: grantbuster Date: Fri, 3 Mar 2023 13:48:42 -0700 Subject: [PATCH 3/4] incremented version for new bias correction feature --- rex/version.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/rex/version.py b/rex/version.py index 6a1a1dad..8ef1da0f 100644 --- a/rex/version.py +++ b/rex/version.py @@ -1,3 +1,3 @@ """rex Version number""" -__version__ = "0.2.78" +__version__ = "0.2.79" From a2a483082a4755196881b9a34079906734457ff8 Mon Sep 17 00:00:00 2001 From: grantbuster Date: Tue, 7 Mar 2023 08:56:50 -0700 Subject: [PATCH 4/4] PR edits --- rex/sam_resource.py | 29 +++++++++++++++-------------- 1 file changed, 15 insertions(+), 14 deletions(-) diff --git a/rex/sam_resource.py b/rex/sam_resource.py index f2e2c6d8..990e49a3 100644 --- a/rex/sam_resource.py +++ b/rex/sam_resource.py @@ -707,12 +707,12 @@ def bias_correct(self, bc_df): Parameters ---------- - bc_df : None | pd.DataFrame - None if not provided or extracted DataFrame with wind or solar - resource bias correction table. This has columns: gid, adder, - scalar. If either adder or scalar is not present, scalar defaults - to 1 and adder to 0. Only windspeed, GHI, and DNI are corrected. - GHI and DNI are corrected with the same correction factors. + bc_df : pd.DataFrame + DataFrame with wind or solar resource bias correction table. This + has columns: gid, adder, scalar. If either adder or scalar is not + present, scalar defaults to 1 and adder to 0. Only windspeed or + GHI+DNI are corrected depending on the technology. GHI and DNI are + corrected with the same correction factors. """ if 'adder' not in bc_df: @@ -743,14 +743,7 @@ def bias_correct(self, bc_df): scalar = np.expand_dims(bc_df.loc[site_arr[isin], 'scalar'].values, 0) adder = np.expand_dims(bc_df.loc[site_arr[isin], 'adder'].values, 0) - if 'windspeed' in self._res_arrays: - logger.debug('Bias correcting windspeeds.') - arr = self._res_arrays['windspeed'] - arr[:, isin] = arr[:, isin] * scalar + adder - arr = np.maximum(arr, 0) - self._res_arrays['windspeed'] = arr - - elif 'ghi' in self._res_arrays: + if 'ghi' in self._res_arrays and 'dni' in self._res_arrays: ghi = self._res_arrays['ghi'] dni = self._res_arrays['dni'] dhi = self._res_arrays['dhi'] @@ -777,7 +770,15 @@ def bias_correct(self, bc_df): self._res_arrays['dni'] = dni self._res_arrays['dhi'] = dhi + elif 'windspeed' in self._res_arrays: + logger.debug('Bias correcting windspeeds.') + arr = self._res_arrays['windspeed'] + arr[:, isin] = arr[:, isin] * scalar + adder + arr = np.maximum(arr, 0) + self._res_arrays['windspeed'] = arr + if self._mean_arrays is not None: + # pylint: disable=consider-iterating-dictionary for var in self._mean_arrays.keys(): self._mean_arrays[var] = self._res_arrays[var].mean(axis=0)