diff --git a/docs/sphinx/source/whatsnew/v0.7.0.rst b/docs/sphinx/source/whatsnew/v0.7.0.rst index 505bd00825..dd2e580413 100644 --- a/docs/sphinx/source/whatsnew/v0.7.0.rst +++ b/docs/sphinx/source/whatsnew/v0.7.0.rst @@ -20,6 +20,10 @@ Bug fixes * Fix handling of keyword arguments in `forecasts.get_processed_data`. (:issue:`745`) * Fix output as Series feature in :py:func:`pvlib.pvsystem.ashraeiam`. +* Fix rounding issue in `clearsky._linearly_scale`, a function that converts + longitude or latitude degree to an index number in a Linke turbidity lookup + table. Also rename the function to `clearsky._degrees_to_index`. + (:issue:`754`) Testing ~~~~~~~ diff --git a/pvlib/clearsky.py b/pvlib/clearsky.py index 0aa80129d9..073ee7844a 100644 --- a/pvlib/clearsky.py +++ b/pvlib/clearsky.py @@ -149,9 +149,9 @@ def lookup_linke_turbidity(time, latitude, longitude, filepath=None, ---------- time : pandas.DatetimeIndex - latitude : float + latitude : float or int - longitude : float + longitude : float or int filepath : None or string, default None The path to the ``.h5`` file. @@ -193,22 +193,12 @@ def lookup_linke_turbidity(time, latitude, longitude, filepath=None, pvlib_path = os.path.dirname(os.path.abspath(__file__)) filepath = os.path.join(pvlib_path, 'data', 'LinkeTurbidities.h5') - latitude_index = ( - np.around(_linearly_scale(latitude, 90, -90, 0, 2160)) - .astype(np.int64)) - longitude_index = ( - np.around(_linearly_scale(longitude, -180, 180, 0, 4320)) - .astype(np.int64)) + latitude_index = _degrees_to_index(latitude, coordinate='latitude') + longitude_index = _degrees_to_index(longitude, coordinate='longitude') - lt_h5_file = tables.open_file(filepath) - try: + with tables.open_file(filepath) as lt_h5_file: lts = lt_h5_file.root.LinkeTurbidity[latitude_index, longitude_index, :] - except IndexError: - raise IndexError('Latitude should be between 90 and -90, ' - 'longitude between -180 and 180.') - finally: - lt_h5_file.close() if interp_turbidity: linke_turbidity = _interpolate_turbidity(lts, time) @@ -299,28 +289,65 @@ def _calendar_month_middles(year): return middles -def _linearly_scale(inputmatrix, inputmin, inputmax, outputmin, outputmax): - """linearly scale input to output, used by Linke turbidity lookup""" +def _degrees_to_index(degrees, coordinate): + """Transform input degrees to an output index integer. The Linke + turbidity lookup tables have three dimensions, latitude, longitude, and + month. Specify a degree value and either 'latitude' or 'longitude' to get + the appropriate index number for the first two of these index numbers. + + Parameters + ---------- + degrees : float or int + Degrees of either latitude or longitude. + coordinate : string + Specify whether degrees arg is latitude or longitude. Must be set to + either 'latitude' or 'longitude' or an error will be raised. + + Returns + ------- + index : np.int16 + The latitude or longitude index number to use when looking up values + in the Linke turbidity lookup table. + """ + # Assign inputmin, inputmax, and outputmax based on degree type. + if coordinate == 'latitude': + inputmin = 90 + inputmax = -90 + outputmax = 2160 + elif coordinate == 'longitude': + inputmin = -180 + inputmax = 180 + outputmax = 4320 + else: + raise IndexError("coordinate must be 'latitude' or 'longitude'.") + inputrange = inputmax - inputmin - outputrange = outputmax - outputmin - delta = outputrange/inputrange # number of indices per input unit - inputmin = inputmin + 1.0 / delta / 2.0 # shift to center of index - outputmax = outputmax - 1 # shift index to zero indexing - outputmatrix = (inputmatrix - inputmin) * delta + outputmin + scale = outputmax/inputrange # number of indices per degree + center = inputmin + 1 / scale / 2 # shift to center of index + outputmax -= 1 # shift index to zero indexing + index = (degrees - center) * scale err = IndexError('Input, %g, is out of range (%g, %g).' % - (inputmatrix, inputmax - inputrange, inputmax)) - # round down if input is within half an index or else raise index error - if outputmatrix > outputmax: - if np.around(outputmatrix - outputmax, 1) <= 0.5: - outputmatrix = outputmax + (degrees, inputmin, inputmax)) + + # If the index is still out of bounds after rounding, raise an error. + # 0.500001 is used in comparisons instead of 0.5 to allow for a small + # margin of error which can occur when dealing with floating point numbers. + if index > outputmax: + if index - outputmax <= 0.500001: + index = outputmax else: raise err - elif outputmatrix < outputmin: - if np.around(outputmin - outputmatrix, 1) <= 0.5: - outputmatrix = outputmin + elif index < 0: + if -index <= 0.500001: + index = 0 else: raise err - return outputmatrix + # If the index wasn't set to outputmax or 0, round it and cast it as an + # integer so it can be used in integer-based indexing. + else: + index = int(np.around(index)) + + return index def haurwitz(apparent_zenith): diff --git a/pvlib/test/test_clearsky.py b/pvlib/test/test_clearsky.py index ee291345f9..44d2e586a7 100644 --- a/pvlib/test/test_clearsky.py +++ b/pvlib/test/test_clearsky.py @@ -519,6 +519,13 @@ def monthly_lt_nointerp(lat, lon, time=months): monthly_lt_nointerp(38.2, -181) # exceeds min longitude +def test_degrees_to_index_1(): + """Test that _degrees_to_index raises an error when something other than + 'latitude' or 'longitude' is passed.""" + with pytest.raises(IndexError): # invalid value for coordinate argument + clearsky._degrees_to_index(degrees=22.0, coordinate='width') + + @pytest.fixture def detect_clearsky_data(): test_dir = os.path.dirname(os.path.abspath(