Skip to content

fix read_tmy3 with year coerced not monotonic, breaks soiling #910

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 15 commits into from
Feb 29, 2020
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
35 changes: 19 additions & 16 deletions docs/examples/plot_greensboro_kimber_soiling.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,20 +2,20 @@
Kimber Soiling Model
====================

Examples of soiling using the Kimber model [1]_.

References
----------
.. [1] "The Effect of Soiling on Large Grid-Connected Photovoltaic Systems
in California and the Southwest Region of the United States," Adrianne
Kimber, et al., IEEE 4th World Conference on Photovoltaic Energy
Conference, 2006, :doi:`10.1109/WCPEC.2006.279690`
Examples of soiling using the Kimber model.
"""

# %%
# This example shows basic usage of pvlib's Kimber Soiling model with
# This example shows basic usage of pvlib's Kimber Soiling model [1]_ with
# :py:meth:`pvlib.losses.soiling_kimber`.
#
# References
# ----------
# .. [1] "The Effect of Soiling on Large Grid-Connected Photovoltaic Systems
# in California and the Southwest Region of the United States," Adrianne
# Kimber, et al., IEEE 4th World Conference on Photovoltaic Energy
# Conference, 2006, :doi:`10.1109/WCPEC.2006.279690`
#
# The Kimber Soiling model assumes that soiling builds up at a constant rate
# until cleaned either manually or by rain. The rain must reach a threshold to
# clean the panels. When rains exceeds the threshold, it's assumed the earth is
Expand All @@ -30,19 +30,22 @@
# step.

from datetime import datetime
import pathlib
from matplotlib import pyplot as plt
from pvlib.iotools import read_tmy3
from pvlib.losses import soiling_kimber
from pvlib.tests.conftest import DATA_DIR
import pvlib

# get full path to the data directory
DATA_DIR = pathlib.Path(pvlib.__file__).parent / 'data'

# get TMY3 data with rain
greensboro = read_tmy3(DATA_DIR / '723170TYA.CSV', coerce_year=1990)
# NOTE: can't use Sand Point, AK b/c Lprecipdepth is -9900, ie: missing
greensboro_rain = greensboro[0].Lprecipdepth
# calculate soiling with no wash dates
greensboro, _ = read_tmy3(DATA_DIR / '723170TYA.CSV', coerce_year=1990)
# get the rain data
greensboro_rain = greensboro.Lprecipdepth
# calculate soiling with no wash dates and cleaning threshold of 25-mm of rain
THRESHOLD = 25.0
soiling_no_wash = soiling_kimber(
greensboro_rain, cleaning_threshold=THRESHOLD, istmy=True)
soiling_no_wash = soiling_kimber(greensboro_rain, cleaning_threshold=THRESHOLD)
soiling_no_wash.name = 'soiling'
# daily rain totals
daily_rain = greensboro_rain.iloc[:-1].resample('D').sum()
Expand Down
21 changes: 13 additions & 8 deletions docs/sphinx/source/whatsnew/v0.7.2.rst
Original file line number Diff line number Diff line change
Expand Up @@ -7,42 +7,47 @@ API Changes
~~~~~~~~~~~
* :py:class:`pvlib.forecast.ForecastModel` now requires ``start`` and ``end``
arguments to be tz-localized. (:issue:`877`, :pull:`879`)
* :py:func:`pvlib.iotools.read_tmy3` when coerced to a single year now returns
indices that are monotonically increasing. Therefore, the last index will be
January 1, 00:00 of the *next* year. (:pull:`910`)

Enhancements
~~~~~~~~~~~~
* TMY3 dataframe returned by :py:func:`~pvlib.iotools.read_tmy3` now contains
the original ``Date (MM/DD/YYYY)`` and ``Time (HH:MM)`` columns that the
indices were parsed from (:pull:`866`)
indices were parsed from. (:pull:`866`)
* Add :py:func:`~pvlib.pvsystem.PVSystem.faiman` and added
``temperature_model='faiman'`` option to :py:class:`~pvlib.modelchain.ModelChain`
(:pull:`897`) (:issue:`836`).
* Add Kimber soiling model :py:func:`pvlib.losses.soiling_kimber` (:pull:`860`)
* Add Kimber soiling model :py:func:`pvlib.losses.soiling_kimber`. (:pull:`860`)

Bug fixes
~~~~~~~~~
* Fix :py:func:`~pvlib.iotools.read_tmy3` parsing when February contains
a leap year (:pull:`866`)
a leap year. (:pull:`866`)
* Implement NREL Developer Network API key for consistent success with API
calls in :py:mod:`pvlib.tests.iotools.test_psm3` (:pull:`873`)
calls in :py:mod:`pvlib.tests.iotools.test_psm3`. (:pull:`873`)
* Fix issue with :py:class:`pvlib.location.Location` creation when
passing ``tz=datetime.timezone.utc`` (:pull:`879`)
passing ``tz=datetime.timezone.utc``. (:pull:`879`)
* Fix documentation homepage title to "pvlib python" based on first heading on
the page. (:pull:`890`) (:issue:`888`)
* Implement `pytest-remotedata <https://github.com/astropy/pytest-remotedata>`_
to increase test suite speed. Requires ``--remote-data`` pytest flag to
execute data retrieval tests over a network.(:issue:`882`)(:pull:`896`)
execute data retrieval tests over a network. (:issue:`882`)(:pull:`896`)
* Fix missing
`0.7.0 what's new <https://pvlib-python.readthedocs.io/en/stable/whatsnew.html#v0-7-0-december-18-2019>`_
entries about changes to ``PVSystem.pvwatts_ac``. Delete unreleased
0.6.4 what's new file. (:issue:`898`)
* Compatibility with cftime 1.1. (:issue:`895`)
* Add Python3.8 to Azure Pipelines CI (:issue:`903`)(:pull:`904`)
* Add documentation build test to Azure Pipelines CI (:pull:`909`)
* Add Python3.8 to Azure Pipelines CI. (:issue:`903`)(:pull:`904`)
* Add documentation build test to Azure Pipelines CI. (:pull:`909`)
* Minor implemention changes to avoid runtime and deprecation warnings in
:py:func:`~pvlib.clearsky.detect_clearsky`,
:py:func:`~pvlib.iam.martin_ruiz_diffuse`,
:py:func:`~pvlib.losses.soiling_hsu`,
and various test functions.
* Fix :py:func:`~pvlib.iotools.read_tmy3` so that when coerced to a single year
the TMY3 index will be monotonically increasing. (:pull:`910`)

Testing
~~~~~~~
Expand Down
3 changes: 1 addition & 2 deletions pvlib/iotools/__init__.py
Original file line number Diff line number Diff line change
@@ -1,5 +1,4 @@
from pvlib.iotools.tmy import read_tmy2 # noqa: F401
from pvlib.iotools.tmy import read_tmy3 # noqa: F401
from pvlib.iotools.tmy import read_tmy2, read_tmy3 # noqa: F401
from pvlib.iotools.epw import read_epw, parse_epw # noqa: F401
from pvlib.iotools.srml import read_srml # noqa: F401
from pvlib.iotools.srml import read_srml_month_from_solardat # noqa: F401
Expand Down
27 changes: 15 additions & 12 deletions pvlib/iotools/tmy.py
Original file line number Diff line number Diff line change
Expand Up @@ -28,10 +28,12 @@ def read_tmy3(filename=None, coerce_year=None, recolumn=True):
a relative file path, absolute file path, or url.

coerce_year : None or int, default None
If supplied, the year of the data will be set to this value.
If supplied, the year of the index will be set to `coerce_year`, except
for the last index value which will be set to the *next* year so that
the index increases monotonically.

recolumn : bool, default True
If True, apply standard names to TMY3 columns. Typically this
If ``True``, apply standard names to TMY3 columns. Typically this
results in stripping the units from the column name.

Returns
Expand All @@ -49,7 +51,6 @@ def read_tmy3(filename=None, coerce_year=None, recolumn=True):

Notes
-----

The returned structures have the following fields.

=============== ====== ===================
Expand Down Expand Up @@ -139,15 +140,16 @@ def read_tmy3(filename=None, coerce_year=None, recolumn=True):
TMYData.PresWthUncertainty Present weather code uncertainty, see [2]_.
============================= ======================================================================================================================================================

.. warning:: TMY3 irradiance data corresponds to the previous hour, so the
first hour is 1AM, corresponding to the net irradiance from midnite to
1AM, and the last hour is midnite of the *next* year, unless the year
has been coerced. EG: if TMY3 was 1988-12-31 24:00:00 this becomes
1989-01-01 00:00:00
.. warning:: TMY3 irradiance data corresponds to the *previous* hour, so
the first index is 1AM, corresponding to the irradiance from midnight
to 1AM, and the last index is midnight of the *next* year. For example,
if the last index in the TMY3 file was 1988-12-31 24:00:00 this becomes
1989-01-01 00:00:00 after calling :func:`~pvlib.iotools.read_tmy3`.

.. warning:: When coercing the year, the last index in the dataframe will
be the first hour of the same year, EG: if TMY3 was 1988-12-31 24:00:00
and year is coerced to 1990 this becomes 1990-01-01
become midnight of the *next* year. For example, if the last index in
the TMY3 was 1988-12-31 24:00:00, and year is coerced to 1990 then this
becomes 1991-01-01 00:00:00.

References
----------
Expand Down Expand Up @@ -214,11 +216,12 @@ def read_tmy3(filename=None, coerce_year=None, recolumn=True):
data_ymd[leapday] += datetime.timedelta(days=1)
# shifted_hour is a pd.Series, so use pd.to_timedelta to get a pd.Series of
# timedeltas
if coerce_year is not None:
data_ymd = data_ymd.map(lambda dt: dt.replace(year=coerce_year))
data_ymd.iloc[-1] = data_ymd.iloc[-1].replace(year=coerce_year+1)
# NOTE: as of pvlib-0.6.3, min req is pandas-0.18.1, so pd.to_timedelta
# unit must be in (D,h,m,s,ms,us,ns), but pandas>=0.24 allows unit='hour'
data.index = data_ymd + pd.to_timedelta(shifted_hour, unit='h')
if coerce_year is not None:
data.index = data.index.map(lambda dt: dt.replace(year=coerce_year))

if recolumn:
data = _recolumn(data) # rename to standard column names
Expand Down
46 changes: 21 additions & 25 deletions pvlib/losses.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,11 +11,14 @@


def soiling_hsu(rainfall, cleaning_threshold, tilt, pm2_5, pm10,
depo_veloc={'2_5': 0.004, '10': 0.0009},
rain_accum_period=pd.Timedelta('1h')):
depo_veloc=None, rain_accum_period=pd.Timedelta('1h')):
"""
Calculates soiling ratio given particulate and rain data using the model
from Humboldt State University [1]_.
from Humboldt State University (HSU).

The HSU soiling model [1]_ returns the soiling ratio, a value between zero
and one which is equivalent to (1 - transmission loss). Therefore a soiling
ratio of 1.0 is equivalent to zero transmission loss.

Parameters
----------
Expand Down Expand Up @@ -65,6 +68,10 @@ def soiling_hsu(rainfall, cleaning_threshold, tilt, pm2_5, pm10,
except ImportError:
raise ImportError("The soiling_hsu function requires scipy.")

# never use mutable input arguments
if depo_veloc is None:
depo_veloc = {'2_5': 0.004, '10': 0.0009}

# accumulate rainfall into periods for comparison with threshold
accum_rain = rainfall.rolling(rain_accum_period, closed='right').sum()
# cleaning is True for intervals with rainfall greater than threshold
Expand All @@ -91,12 +98,12 @@ def soiling_hsu(rainfall, cleaning_threshold, tilt, pm2_5, pm10,

def soiling_kimber(rainfall, cleaning_threshold=6, soiling_loss_rate=0.0015,
grace_period=14, max_soiling=0.3, manual_wash_dates=None,
initial_soiling=0, rain_accum_period=24, istmy=False):
initial_soiling=0, rain_accum_period=24):
"""
Calculate soiling ratio with rainfall data and a daily soiling rate using
the Kimber soiling model [1]_.
Calculates fraction of energy lost due to soiling given rainfall data and
daily loss rate using the Kimber model.

Kimber soiling model assumes soiling builds up at a daily rate unless
Kimber soiling model [1]_ assumes soiling builds up at a daily rate unless
the daily rainfall is greater than a threshold. The model also assumes that
if daily rainfall has exceeded the threshold within a grace period, then
the ground is too damp to cause soiling build-up. The model also assumes
Expand Down Expand Up @@ -127,8 +134,6 @@ def soiling_kimber(rainfall, cleaning_threshold=6, soiling_loss_rate=0.0015,
rain_accum_period : int, default 24
Period for accumulating rainfall to check against `cleaning_threshold`.
The Kimber model defines this period as one day. [hours]
istmy : bool, default False
Fix last timestep in TMY so that it is monotonically increasing.

Returns
-------
Expand Down Expand Up @@ -166,23 +171,12 @@ def soiling_kimber(rainfall, cleaning_threshold=6, soiling_loss_rate=0.0015,
# convert grace_period to timedelta
grace_period = datetime.timedelta(days=grace_period)

# get rainfall timezone, timestep as timedelta64, and timestep in int days
rain_tz = rainfall.index.tz
rain_index = rainfall.index.values
timestep_interval = (rain_index[1] - rain_index[0])
# get indices as numpy datetime64, calculate timestep as numpy timedelta64,
# and convert timestep to fraction of days
rain_index_vals = rainfall.index.values
timestep_interval = (rain_index_vals[1] - rain_index_vals[0])
day_fraction = timestep_interval / np.timedelta64(24, 'h')

# if TMY fix to be monotonically increasing by rolling index by 1 interval
# and then adding 1 interval, while the values stay the same
if istmy:
rain_index = np.roll(rain_index, 1) + timestep_interval
# NOTE: numpy datetim64[ns] has no timezone
# convert to datetimeindex at UTC and convert to original timezone
rain_index = pd.DatetimeIndex(rain_index, tz='UTC').tz_convert(rain_tz)
# fixed rainfall timeseries with monotonically increasing index
rainfall = pd.Series(
rainfall.values, index=rain_index, name=rainfall.name)

# accumulate rainfall
accumulated_rainfall = rainfall.rolling(
rain_accum_period, closed='right').sum()
Expand All @@ -191,6 +185,7 @@ def soiling_kimber(rainfall, cleaning_threshold=6, soiling_loss_rate=0.0015,
soiling = np.ones_like(rainfall.values) * soiling_loss_rate * day_fraction
soiling[0] = initial_soiling
soiling = np.cumsum(soiling)
soiling = pd.Series(soiling, index=rainfall.index, name='soiling')

# rainfall events that clean the panels
rain_events = accumulated_rainfall > cleaning_threshold
Expand All @@ -205,8 +200,9 @@ def soiling_kimber(rainfall, cleaning_threshold=6, soiling_loss_rate=0.0015,

# manual wash dates
if manual_wash_dates is not None:
rain_tz = rainfall.index.tz
# convert manual wash dates to datetime index in the timezone of rain
manual_wash_dates = pd.DatetimeIndex(manual_wash_dates, tz=rain_tz)
soiling = pd.Series(soiling, index=rain_index, name='soiling')
cleaning[manual_wash_dates] = soiling[manual_wash_dates]

# remove soiling by foward filling cleaning where NaN
Expand Down
18 changes: 11 additions & 7 deletions pvlib/tests/iotools/test_tmy.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@
import pandas as pd
import pytest
from pvlib.iotools import tmy
from pvlib.iotools import read_tmy3
from conftest import DATA_DIR

# test the API works
Expand Down Expand Up @@ -30,21 +31,23 @@ def test_read_tmy3_recolumn():


def test_read_tmy3_norecolumn():
data, meta = tmy.read_tmy3(TMY3_TESTFILE, recolumn=False)
data, _ = tmy.read_tmy3(TMY3_TESTFILE, recolumn=False)
assert 'GHI source' in data.columns


def test_read_tmy3_coerce_year():
coerce_year = 1987
data, meta = tmy.read_tmy3(TMY3_TESTFILE, coerce_year=coerce_year)
assert (data.index.year == 1987).all()
data, _ = tmy.read_tmy3(TMY3_TESTFILE, coerce_year=coerce_year)
assert (data.index[:-1].year == 1987).all()
assert data.index[-1].year == 1988


def test_read_tmy3_no_coerce_year():
coerce_year = None
data, meta = tmy.read_tmy3(TMY3_TESTFILE, coerce_year=coerce_year)
data, _ = tmy.read_tmy3(TMY3_TESTFILE, coerce_year=coerce_year)
assert 1997 and 1999 in data.index.year

assert data.index[-2] == pd.Timestamp('1998-12-31 23:00:00-09:00')
assert data.index[-1] == pd.Timestamp('1999-01-01 00:00:00-09:00')

def test_read_tmy2():
tmy.read_tmy2(TMY2_TESTFILE)
Expand All @@ -70,9 +73,10 @@ def test_gh865_read_tmy3_feb_leapyear_hr24():
# now check if it parses correctly when we try to coerce the year
data, _ = read_tmy3(TMY3_FEB_LEAPYEAR, coerce_year=1990)
# if get's here w/o an error, then gh865 is fixed, but let's check anyway
assert all(data.index.year == 1990)
assert all(data.index[:-1].year == 1990)
assert data.index[-1].year == 1991
# let's do a quick sanity check, are the indices monotonically increasing?
assert all(np.diff(data.index[:-1].astype(int)) == 3600000000000)
assert all(np.diff(data.index.astype(int)) == 3600000000000)
# according to the TMY3 manual, each record corresponds to the previous
# hour so check that the 1st hour is 1AM and the last hour is midnite
assert data.index[0].hour == 1
Expand Down
Loading