[DEV] address pandas 3.1.0 compatibility issues

benchmarks/benchmarks/solarposition.py

@@ -51,7 +51,7 @@ def time_sun_rise_set_transit_ephem(self, ndays):
             self.times_daily, self.lat, self.lon)
 
     def time_sun_rise_set_transit_geometric_full_comparison(self, ndays):
-        dayofyear = self.times_daily.dayofyear
+        dayofyear = self.times_daily.day_of_year
         declination = solarposition.declination_spencer71(dayofyear)
         equation_of_time = solarposition.equation_of_time_spencer71(dayofyear)
         solarposition.sun_rise_set_transit_geometric(

docs/examples/solar-position/plot_sunpath_diagrams.py

@@ -32,13 +32,13 @@
 ax = plt.subplot(1, 1, 1, projection='polar')
 # draw the analemma loops
 points = ax.scatter(np.radians(solpos.azimuth), solpos.apparent_zenith,
-                    s=2, label=None, c=solpos.index.dayofyear,
+                    s=2, label=None, c=solpos.index.day_of_year,
                     cmap='twilight_shifted_r')
 # add and format colorbar
 cbar = ax.figure.colorbar(points)
 times_ticks = pd.date_range('2019-01-01', '2020-01-01', freq='MS', tz=tz)
-cbar.set_ticks(ticks=times_ticks.dayofyear, labels=[], minor=False)
-cbar.set_ticks(ticks=times_ticks.dayofyear+15,
+cbar.set_ticks(ticks=times_ticks.day_of_year, labels=[], minor=False)
+cbar.set_ticks(ticks=times_ticks.day_of_year+15,
                labels=times_ticks.strftime('%b'),
                minor=True)
 cbar.ax.tick_params(which='minor', width=0)
@@ -120,13 +120,13 @@
 
 fig, ax = plt.subplots()
 points = ax.scatter(solpos.azimuth, solpos.apparent_elevation, s=2,
-                    c=solpos.index.dayofyear, label=None,
+                    c=solpos.index.day_of_year, label=None,
                     cmap='twilight_shifted_r')
 # add and format colorbar
 cbar = fig.colorbar(points)
 times_ticks = pd.date_range('2019-01-01', '2020-01-01', freq='MS', tz=tz)
-cbar.set_ticks(ticks=times_ticks.dayofyear, labels=[], minor=False)
-cbar.set_ticks(ticks=times_ticks.dayofyear+15,
+cbar.set_ticks(ticks=times_ticks.day_of_year, labels=[], minor=False)
+cbar.set_ticks(ticks=times_ticks.day_of_year+15,
                labels=times_ticks.strftime('%b'),
                minor=True)
 cbar.ax.tick_params(which='minor', width=0)

docs/sphinx/source/whatsnew/v0.16.2.rst

@@ -24,6 +24,7 @@ Bug fixes
 
 Enhancements
 ~~~~~~~~~~~~
+* Add ``pandas==3.1.0`` support. (:pull:`2924`)
 * Add ``dataset`` parameter to :py:func:`~pvlib.iotools.get_era5`
   allowing for choosing between ERA5 and ERA5-Land datasets.
   (:pull:`2780`)

pvlib/clearsky.py

@@ -245,7 +245,7 @@ def _interpolate_turbidity(lts, time):
 
     isleap = time_utc.is_leap_year
 
-    dayofyear = time_utc.dayofyear
+    dayofyear = time_utc.day_of_year
     days_leap = _calendar_month_middles(2016)
     days_no_leap = _calendar_month_middles(2015)
 

pvlib/iotools/pvgis.py

@@ -766,7 +766,7 @@ def get_pvgis_horizon(latitude, longitude, url=URL, **kwargs):
     data.columns = ['horizon_azimuth', 'horizon_elevation']
     # Convert azimuth to pvlib convention (north=0, south=180)
     data['horizon_azimuth'] += 180
-    data.set_index('horizon_azimuth', inplace=True)
+    data = data.set_index('horizon_azimuth')
     data = data['horizon_elevation']  # convert to pd.Series
     data = data[data.index < 360]  # remove duplicate north point (0 and 360)
     return data, metadata

pvlib/iotools/solargis.py

@@ -204,9 +204,8 @@ def get_solargis(latitude, longitude, start, end, variables, api_key,
     if map_variables:
         for variable in VARIABLE_MAP:
             if variable.solargis_name in data.columns:
-                data.rename(
-                    columns={variable.solargis_name: variable.pvlib_name},
-                    inplace=True
+                data = data.rename(
+                    columns={variable.solargis_name: variable.pvlib_name}
                 )
                 data[variable.pvlib_name] = data[
                     variable.pvlib_name].apply(variable.conversion)

pvlib/iotools/solcast.py

@@ -421,9 +421,8 @@ def _solcast2pvlib(data):
     # rename and convert variables
     for variable in VARIABLE_MAP:
         if variable.solcast_name in data.columns:
-            data.rename(
+            data = data.rename(
                 columns={variable.solcast_name: variable.pvlib_name},
-                inplace=True
             )
             data[variable.pvlib_name] = data[
                 variable.pvlib_name].apply(variable.conversion)
@@ -474,7 +473,7 @@ def _get_solcast(
 
     if response.status_code == 200:
         j = response.json()
-        df = pd.DataFrame.from_dict(j[list(j.keys())[0]])
+        df = pd.DataFrame(j[list(j.keys())[0]])
         if map_variables:
             return _solcast2pvlib(df)
         else:

pvlib/iotools/solrad.py

@@ -198,7 +198,7 @@
     dates = pd.date_range(start.floor('D'), end, freq='D')
     station = station.lower()
     filenames = [
-        f"{station}/{d.year}/{station}{d.strftime('%y')}{d.dayofyear:03}.dat"
+        f"{station}/{d.year}/{station}{d.strftime('%y')}{d.day_of_year:03}.dat"
         for d in dates
     ]
 

pvlib/iotools/surfrad.py

@@ -158,7 +158,7 @@ def read_surfrad(filename, map_variables=True):
     data = data.where(~missing, np.nan)
 
     if map_variables:
-        data.rename(columns=VARIABLE_MAP, inplace=True)
+        data = data.rename(columns=VARIABLE_MAP)
     return data, metadata
 
 
@@ -250,7 +250,7 @@ def get_surfrad(station, start, end, map_variables=True,
     station = station.lower()
 
     filenames = [
-        f"{station}/{d.year}/{station}{d.strftime('%y')}{d.dayofyear:03}.dat"  # noqa: E231,E501
+        f"{station}/{d.year}/{station}{d.strftime('%y')}{d.day_of_year:03}.dat"  # noqa: E231,E501
         for d in dates
     ]
 

pvlib/iotools/tmy.py

@@ -222,7 +222,7 @@ def read_tmy3(filename, coerce_year=None, map_variables=True, encoding=None):
     # 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 = data_ymd.transform(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'

pvlib/irradiance.py

@@ -1937,7 +1937,7 @@ def disc(ghi, solar_zenith, datetime_or_doy, pressure=101325,
 
     datetime_or_doy : int, float, array, pd.DatetimeIndex
         Day of year or array of days of year e.g.
-        pd.DatetimeIndex.dayofyear, or pd.DatetimeIndex.
+        pd.DatetimeIndex.day_of_year, or pd.DatetimeIndex.
 
     pressure : numeric or None, default 101325
         Site pressure. See :term:`pressure`. [Pa]
@@ -2765,7 +2765,7 @@ def erbs(ghi, zenith, datetime_or_doy, min_cos_zenith=0.065, max_zenith=87):
         :term:`solar_zenith`. [°]
     datetime_or_doy : int, float, array, pd.DatetimeIndex
         Day of year or array of days of year e.g.
-        pd.DatetimeIndex.dayofyear, or pd.DatetimeIndex.
+        pd.DatetimeIndex.day_of_year, or pd.DatetimeIndex.
     min_cos_zenith : numeric, default 0.065
         Minimum value of cos(zenith) to allow when calculating global
         clearness index `kt`. Equivalent to zenith = 86.273°.
@@ -2865,7 +2865,7 @@ def erbs_driesse(ghi, zenith, datetime_or_doy=None, dni_extra=None,
 
     datetime_or_doy : int, float, array or pd.DatetimeIndex, optional
         Day of year or array of days of year e.g.
-        pd.DatetimeIndex.dayofyear, or pd.DatetimeIndex.
+        pd.DatetimeIndex.day_of_year, or pd.DatetimeIndex.
         Either ``datetime_or_doy`` or ``dni_extra`` must be provided.
 
     dni_extra : numeric, optional
@@ -2992,7 +2992,7 @@ def orgill_hollands(ghi, zenith, datetime_or_doy, dni_extra=None,
 
     datetime_or_doy : int, float, array or pd.DatetimeIndex, optional
         Day of year or array of days of year e.g.
-        pd.DatetimeIndex.dayofyear, or pd.DatetimeIndex.
+        pd.DatetimeIndex.day_of_year, or pd.DatetimeIndex.
         Either ``datetime_or_doy`` or ``dni_extra`` must be provided.
 
     dni_extra : numeric, optional
@@ -3093,7 +3093,7 @@ def boland(ghi, solar_zenith, datetime_or_doy, a_coeff=8.645, b_coeff=0.613,
         :term:`solar_zenith`. [°]
 
     datetime_or_doy : numeric or pd.DatetimeIndex.
-        Day of year or array of days of year e.g. pd.DatetimeIndex.dayofyear,
+        Day of year or array of days of year e.g. pd.DatetimeIndex.day_of_year,
         or pd.DatetimeIndex.
 
     a_coeff : float, default 8.645
@@ -3897,7 +3897,7 @@ def louche(ghi, solar_zenith, datetime_or_doy, max_zenith=90):
 
     datetime_or_doy : numeric, pandas.DatetimeIndex
         Day of year or array of days of year e.g.
-        pd.DatetimeIndex.dayofyear, or pd.DatetimeIndex.
+        pd.DatetimeIndex.day_of_year, or pd.DatetimeIndex.
 
     Returns
     -------

pvlib/ivtools/utils.py

@@ -160,7 +160,7 @@ def rectify_iv_curve(voltage, current, decimals=None):
 
     df = pd.DataFrame(data=np.vstack((voltage, current)).T, columns=['v', 'i'])
     # restrict to first quadrant
-    df.dropna(inplace=True)
+    df = df.dropna()
     df = df[(df['v'] >= 0) & (df['i'] >= 0)]
     # sort pairs on voltage, then current
     df = df.sort_values(by=['v', 'i'], ascending=[True, False])

pvlib/snow.py

@@ -5,7 +5,13 @@
 
 import numpy as np
 import pandas as pd
-from pvlib.tools import sind, cosd, tand
+
+from pvlib.tools import (
+    _pandas_infer_freq_BaseOffset_transition_handler,
+    cosd,
+    sind,
+    tand,
+)
 
 
 def _time_delta_in_hours(times):
@@ -58,7 +64,7 @@ def fully_covered_nrel(snowfall, snow_depth=None, threshold_snowfall=1.,
     hourly_snow_rate = snowfall / timestep
     # if we can infer a time frequency, use first snowfall value
     # otherwise the first snowfall value is ignored
-    freq = pd.infer_freq(snowfall.index)
+    freq = _pandas_infer_freq_BaseOffset_transition_handler(snowfall.index)
     if freq is not None:
         timedelta = pd.tseries.frequencies.to_offset(freq) / pd.Timedelta('1h')
         hourly_snow_rate.iloc[0] = snowfall.iloc[0] / timedelta
@@ -170,7 +176,7 @@ def coverage_nrel(snowfall, poa_irradiance, temp_air, surface_tilt,
     snow_coverage[new_snowfall] = 1.0
     if np.isnan(snow_coverage.iloc[0]):
         snow_coverage.iloc[0] = initial_coverage
-    snow_coverage.ffill(inplace=True)
+    snow_coverage = snow_coverage.ffill()
     snow_coverage -= cumulative_sliding
 
     return snow_coverage.clip(lower=0)

pvlib/soiling.py

@@ -192,9 +192,9 @@ def kimber(rainfall, cleaning_threshold=6, soiling_loss_rate=0.0015,
     # convert grace_period to timedelta
     grace_period = datetime.timedelta(days=grace_period)
 
-    # get indices as numpy datetime64, calculate timestep as numpy timedelta64,
+    # get indices as np.datetime64 (UTC), calculate timestep as np.timedelta64,
     # and convert timestep to fraction of days
-    rain_index_vals = rainfall.index.values
+    rain_index_vals = rainfall.index.tz_convert(None).to_numpy()
     timestep_interval = (rain_index_vals[1] - rain_index_vals[0])
     day_fraction = timestep_interval / np.timedelta64(24, 'h')
 

pvlib/solarposition.py

@@ -772,7 +772,7 @@ def ephemeris(time, latitude, longitude, pressure=101325.0, temperature=12.0):
     time_utc = tools._pandas_to_utc(time)
 
     # strip out the day of the year and calculate the decimal hour
-    DayOfYear = time_utc.dayofyear
+    DayOfYear = time_utc.day_of_year
     DecHours = (time_utc.hour + time_utc.minute/60. + time_utc.second/3600. +
                 time_utc.microsecond/3600.e6)
 

pvlib/tools.py

@@ -9,6 +9,7 @@
 
 import numpy as np
 import pandas as pd
+from packaging.version import Version
 import zoneinfo
 
 
@@ -214,7 +215,7 @@
     -----
     Day of year is determined using UTC, since pandas uses local hour
     """
-    return _pandas_to_utc(pd_object).dayofyear
+    return _pandas_to_utc(pd_object).day_of_year
 
 
 def _pandas_to_utc(pd_object):
@@ -428,10 +429,14 @@
     """ Calculates time interval and samples per window for Reno-style clear
     sky detection functions
     """
-    deltas = np.diff(times.values) / np.timedelta64(1, '60s')
+    deltas = (
+        np.diff(times.tz_localize(None).to_numpy())  # tz-naive np.datetime64
+        / np.timedelta64(1, '60s')
+    )
 
     # determine if we can proceed
-    if times.inferred_freq and len(np.unique(deltas)) == 1:
+    if (_pandas_infer_freq_BaseOffset_transition_handler(times)
+        and len(np.unique(deltas)) == 1):
         sample_interval = times[1] - times[0]
         sample_interval = sample_interval.seconds / 60  # in minutes
         samples_per_window = int(win_length / sample_interval)
@@ -445,7 +450,42 @@
         raise NotImplementedError(message)
 
 
+def _pandas_infer_freq_BaseOffset_transition_handler(times):
+    """
+    Pandas 3.1.0 introduced the following deprecation [^1]:
+        Pandas4Warning: A future version of pandas will return a BaseOffset
+        object instead of a string from infer_freq.
+        Use pd.set_option('future.infer_freq_returns_offset', True) to get the
+        future behavior, or set to False to keep the old behavior and silence
+        this warning.
+        To preserve the string representation, use ``infer_freq(...).freqstr``
+
+    Side-effects must be avoided in user-code.
+
+    pd.infer_freq documentation is at [^2].
+    pd.tseries.offsets.BaseOffset documentation is at [^3].
+
+    Parameters
+    ----------
+    times : DatetimeIndex | TimedeltaIndex | Series | DatetimeLikeArrayMixin
+
+    Returns
+    -------
+    str
+
+    [^1]: pandas PR GH#64679: https://github.com/pandas-dev/pandas/pull/64679
+    [^2]: https://pandas.pydata.org/docs/reference/api/pandas.infer_freq.html
+    [^3]: https://pandas.pydata.org/docs/reference/api/pandas.tseries.offsets.BaseOffset.html
+    """
+    if Version(pd.__version__) < Version("3.1.0rc0"):
+        return pd.infer_freq(times)
+    with pd.option_context({'future.infer_freq_returns_offset': True}):
+        # relies on BaseOffset to work seamlessly with Pandas 4
+        base_offset = pd.infer_freq(times)
+        # base_offset:=None if unable to infer
+        return getattr(base_offset, "freqstr", None)
+
 def _degrees_to_index(degrees, coordinate):
     """Transform input degrees to an output index integer.
     Specify a degree value and either 'latitude' or 'longitude' to get
     the appropriate index number for these two index numbers.

tests/conftest.py

@@ -56,6 +56,10 @@ def _check_pandas_assert_kwargs(kwargs):
 
 def assert_index_equal(left, right, **kwargs):
     kwargs = _check_pandas_assert_kwargs(kwargs)
+    if Version(pd.__version__) >= Version('3.1.0'):
+        # check_freq was added in pandas 3.1.0
+        # tests written prior to it, may have not taken it into account
+        kwargs['check_freq'] = False
     pd.testing.assert_index_equal(left, right, **kwargs)
 
 

tests/iotools/test_era5.py

@@ -8,7 +8,12 @@
 import pvlib
 import requests
 import os
-from tests.conftest import RERUNS, RERUNS_DELAY, requires_ecmwf_credentials
+from tests.conftest import (
+    RERUNS,
+    RERUNS_DELAY,
+    assert_index_equal,
+    requires_ecmwf_credentials,
+)
 
 
 @pytest.fixture
@@ -66,7 +71,7 @@ def test_get_era5_land(params, expected):
     df, meta = pvlib.iotools.get_era5(**params)
     assert np.isclose(meta['longitude'], -80.0)
     assert np.isclose(meta['latitude'], 40.0)
-    pd.testing.assert_index_equal(df.index, expected.index)
+    assert_index_equal(df.index, expected.index)
     assert set(df.columns) == {'ghi', 'temp_air'}
     assert df.notna().all().all()
 

tests/iotools/test_nasa_power.py

@@ -2,7 +2,7 @@
 import pytest
 import pvlib
 from requests.exceptions import HTTPError
-from tests.conftest import RERUNS, RERUNS_DELAY
+from tests.conftest import RERUNS, RERUNS_DELAY, assert_index_equal
 
 
 @pytest.fixture
@@ -39,7 +39,7 @@ def test_get_nasa_power(data_index, ghi_series):
     assert meta['time_standard'] == 'UTC'
     assert meta['title'] == 'NASA/POWER Source Native Resolution Hourly Data'
     # Assert that the index is parsed correctly
-    pd.testing.assert_index_equal(data.index, data_index)
+    assert_index_equal(data.index, data_index)
     # Test one column
     pd.testing.assert_series_equal(data['ALLSKY_SFC_SW_DWN'], ghi_series,
                                    check_freq=False, check_names=False)
@@ -54,7 +54,7 @@ def test_get_nasa_power_pvlib_params_naming(data_index, ghi_series):
                                               end=data_index[-1],
                                               parameters=['ghi'])
     # Assert that the index is parsed correctly
-    pd.testing.assert_index_equal(data.index, data_index)
+    assert_index_equal(data.index, data_index)
     # Test one column
     pd.testing.assert_series_equal(data['ghi'], ghi_series,
                                    check_freq=False)