Repository: pvlib/pvanalytics

Issue:
Function to calculate the Variability Index
I have a function somewhere that calculates the variability index (VI).  I think this package would be a good home for it and I'd be happy to contribute it if there is a place for it.

The features module might be an appropriate place, but a daily float value for VI is returned rather than a boolean.

Discussion:

Patch:
diff --git a/docs/api.rst b/docs/api.rst
index 2ab85a5..39433e9 100644
--- a/docs/api.rst
+++ b/docs/api.rst
@@ -270,4 +270,14 @@ The following functions can be used to calculate system performance metrics.
 .. autosummary::
    :toctree: generated/
 
-   metrics.performance_ratio_nrel
\ No newline at end of file
+   metrics.performance_ratio_nrel
+
+Variability
+-----------
+
+Functions to calculate variability statistics.
+
+.. autosummary::
+   :toctree: generated/
+
+   metrics.variability_index
diff --git a/docs/whatsnew/0.1.1.rst b/docs/whatsnew/0.1.1.rst
new file mode 100644
index 0000000..e994726
--- /dev/null
+++ b/docs/whatsnew/0.1.1.rst
@@ -0,0 +1,22 @@
+.. _whatsnew_011:
+
+0.1.1 (MONTH DAY, YEAR)
+-------------------------
+
+API Changes
+~~~~~~~~~~~
+
+Enhancements
+~~~~~~~~~~~~
+
+* Quantification of irradiance variability with
+  :py:func:`pvanalytics.metrics.variability_index`.
+
+
+Bug Fixes
+~~~~~~~~~
+
+Contributors
+~~~~~~~~~~~~
+
+* Kevin Anderson (:ghuser:`kanderso-nrel`)
diff --git a/pvanalytics/metrics.py b/pvanalytics/metrics.py
index eb7601c..bac36a2 100644
--- a/pvanalytics/metrics.py
+++ b/pvanalytics/metrics.py
@@ -65,3 +65,43 @@ def performance_ratio_nrel(poa_global, temp_air, wind_speed, pac, pdc0,
     performance_ratio = pac.sum() / pdc.sum()
 
     return performance_ratio
+
+
+def _calc_pathlength(signal, freq):
+    # utility function to calculate the arc length of a time series.
+    # used when calculating the variability index.
+    dt = signal.index.to_series(keep_tz=True).diff().dt.total_seconds()/60
+    dy = signal.diff()
+    d = (dy**2 + dt**2)**0.5
+    if freq is not None:
+        return d.resample(freq).sum()
+    return d.sum()
+
+
+def variability_index(measured, clearsky, freq=None):
+    """
+    Calculate the variability index.
+
+    Parameters
+    ----------
+    measured : Series
+        Time series of measured GHI. [W/m2]
+    clearsky : Series
+        Time series of the expected clearsky GHI. [W/m2]
+    freq : pandas datetime offset, optional
+        Aggregation period (e.g. 'D' for daily).  If not specified,
+        the variability index for the entire time series will be returned.
+
+    Returns
+    -------
+    vi : Series or float
+        The calculated variability index
+
+    References
+    ----------
+    .. [1] Stein, Joshua, Hansen, Clifford, and Reno, Matthew J. THE
+       VARIABILITY INDEX: A NEW AND NOVEL METRIC FOR QUANTIFYING IRRADIANCE
+       AND PV OUTPUT VARIABILITY. SAND2012-2088C, World Renewable Energy Forum,
+       2012.
+    """
+    return _calc_pathlength(measured, freq) / _calc_pathlength(clearsky, freq)


Test patch:
diff --git a/.github/workflows/lint-and-test.yml b/.github/workflows/lint-and-test.yml
index 5b1bcb2..47a0a50 100644
--- a/.github/workflows/lint-and-test.yml
+++ b/.github/workflows/lint-and-test.yml
@@ -39,6 +39,7 @@ jobs:
         GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }}
         COVERALLS_FLAG_NAME: ${{ matrix.os }}-${{ matrix.python-version }}-${{ matrix.requirements }}
         COVERALLS_PARALLEL: true
+        COVERALLS_SERVICE_NAME: github
 
   coveralls:
 
diff --git a/pvanalytics/tests/features/test_clearsky.py b/pvanalytics/tests/features/test_clearsky.py
index 8bbfab5..f32953d 100644
--- a/pvanalytics/tests/features/test_clearsky.py
+++ b/pvanalytics/tests/features/test_clearsky.py
@@ -4,6 +4,7 @@ import pandas as pd
 from pvanalytics.features import clearsky
 
 
+@pytest.mark.xfail(strict=True)
 @pytest.mark.filterwarnings("ignore:Support for multi-dimensional indexing")
 def test_reno_identical(quadratic):
     """Identical clearsky and measured irradiance all True"""
@@ -13,6 +14,7 @@ def test_reno_identical(quadratic):
     assert clearsky.reno(quadratic, quadratic).all()
 
 
+@pytest.mark.xfail(strict=True)
 @pytest.mark.filterwarnings("ignore:Support for multi-dimensional indexing")
 @pytest.mark.filterwarnings("ignore:invalid value encountered in")
 def test_reno_begining_end(quadratic):
diff --git a/pvanalytics/tests/test_metrics.py b/pvanalytics/tests/test_metrics.py
index 5438233..36935b3 100644
--- a/pvanalytics/tests/test_metrics.py
+++ b/pvanalytics/tests/test_metrics.py
@@ -1,4 +1,5 @@
 import numpy as np
+import pandas as pd
 from pvanalytics import metrics
 import pytest
 
@@ -23,3 +24,42 @@ def test_performance_ratio_nrel():
                                                        wind_speed, pac, pdc0,
                                                        a, b, deltaT, gamma_pdc)
     assert performance_ratio == pytest.approx(expected)
+
+
+@pytest.fixture
+def variability_inputs():
+    times = pd.date_range('2019-01-01', freq='1min', periods=60*24)
+    clear = pd.Series(100.0, index=times)
+    # alternating sawtooth:
+    jagged = pd.Series([100.0, 101.0] * (len(times)//2), index=times)
+    return pd.DataFrame({
+        'clear': clear,
+        'jagged': jagged
+    })
+
+
+def test_variability_index(variability_inputs):
+    # default freq parameter
+    clear = variability_inputs['clear']
+    jagged = variability_inputs['jagged']
+
+    clear_clear = metrics.variability_index(clear, clear)
+    assert clear_clear == pytest.approx(1.0)
+
+    jagged_clear = metrics.variability_index(jagged, clear)
+    assert jagged_clear == pytest.approx(2**0.5)
+
+
+def test_variability_index_freq(variability_inputs):
+    # custom freq parameter
+    clear = variability_inputs['clear']
+    jagged = variability_inputs['jagged']
+    times = clear.resample('h').mean().index
+
+    expected_clear_clear = pd.Series(1.0, index=times)
+    actual_clear_clear = metrics.variability_index(clear, clear, freq='h')
+    pd.testing.assert_series_equal(actual_clear_clear, expected_clear_clear)
+
+    expected_jagged_clear = pd.Series(2**0.5, index=times)
+    actual_jagged_clear = metrics.variability_index(jagged, clear, freq='h')
+    pd.testing.assert_series_equal(actual_jagged_clear, expected_jagged_clear)