# Source code for pytesmo.utils

# Copyright (c) 2015,Vienna University of Technology,
# Department of Geodesy and Geoinformation

# Redistribution and use in source and binary forms, with or without
# modification, are permitted provided that the following conditions are met:
#   * Redistributions of source code must retain the above copyright
#     notice, this list of conditions and the following disclaimer.
#    * Redistributions in binary form must reproduce the above copyright
#      notice, this list of conditions and the following disclaimer in the
#      documentation and/or other materials provided with the distribution.
#    * Neither the name of the Vienna University of Technology,
#      Department of Geodesy and Geoinformation nor the
#      names of its contributors may be used to endorse or promote products
#      derived from this software without specific prior written permission.

# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
# ARE DISCLAIMED. IN NO EVENT SHALL VIENNA UNIVERSITY OF TECHNOLOGY,
# DEPARTMENT OF GEODESY AND GEOINFORMATION BE LIABLE FOR ANY
# DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
# (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
# LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
# ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
# (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
# SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

'''
Module containing utility functions that do not fit into other modules
'''
import numpy as np
import scipy.interpolate as sc_int
import scipy.optimize as sc_opt
import scipy.special as sc_special

[docs]def ml_percentile(in_data, percentiles):
"""
Calculate percentiles in the way Matlab and IDL do it.

By using interpolation between the lowest an highest rank and the
minimum and maximum outside.

Parameters
----------
in_data: numpy.ndarray
input data
percentiles: numpy.ndarray
percentiles at which to calculate the values

Returns
-------
perc: numpy.ndarray
values of the percentiles
"""

data = np.sort(in_data)
p_rank = 100.0 * (np.arange(data.size) + 0.5) / data.size
perc = np.interp(percentiles, p_rank, data, left=data[0], right=data[-1])
return perc

[docs]def interp_uniq(src):
"""
replace non unique values by their linear interpolated value
This method interpolates iteratively like it is done in IDL.

Parameters
----------
src: numpy.array
array to ensure uniqueness of

Returns
-------
src: numpy.array
interpolated unique values in array of same size as src

"""
size = len(src)
uniq, uniq_ind, counts = np.unique(
src, return_index=True, return_counts=True)

while len(src[uniq_ind]) != size:
# replace non unique percentiles by their linear interpolated value
# This method interpolates iteratively like it is done in IDL
# and might be replaced by a faster method of simple linear
# interpolation
for i in range(len(uniq_ind)):
pos = np.where(src == src[uniq_ind[i]])[0]
if len(pos) > 1:
if pos[0] == 0 and pos[-1] < size - 1:
src[
pos[-1]] = (src[pos[len(pos) - 2]] + src[pos[-1] + 1]) / 2.0
elif pos[-1] == size - 1:
src[pos[0]] = (
src[pos[1]] + src[pos[0] - 1]) / 2.0
else:
src[pos[0]] = (
src[pos[1]] + src[pos[0] - 1]) / 2.0
src[pos[1]] = (
src[pos[0]] + src[pos[1] + 1]) / 2.0
uniq_ind = np.unique(src, return_index=True)[1]

return src

[docs]def unique_percentiles_interpolate(perc_values,
percentiles=[0, 5, 10, 30, 50,
70, 90, 95, 100],
k=1):
"""
Try to ensure that percentile values are unique
and have values for the given percentiles.

If only all the values in perc_values are the same.
The array is unchanged.

Parameters
----------
perc_values: list or numpy.ndarray
calculated values for the given percentiles
percentiles: list or numpy.ndarray
Percentiles to use for CDF matching
k: int
Degree of spline interpolation to use for
filling duplicate percentile values

Returns
-------
uniq_perc_values: numpy.ndarray
Unique percentile values generated through linear
interpolation over removed duplicate percentile values
"""
uniq_ind = np.unique(perc_values, return_index=True)[1]
if len(uniq_ind) == 1:
uniq_ind = np.repeat(uniq_ind, 2)
uniq_ind[-1] = len(percentiles) - 1
uniq_perc_values = perc_values[uniq_ind]

inter = sc_int.InterpolatedUnivariateSpline(
np.array(percentiles)[uniq_ind],
uniq_perc_values,
k=k, ext=0,
check_finite=True)
uniq_perc_values = inter(percentiles)
return uniq_perc_values

[docs]def unique_percentiles_beta(perc_values,
percentiles):
"""
Compute unique percentile values
by fitting the CDF of a beta distribution to the
percentiles.

Parameters
----------
perc_values: list or numpy.ndarray
calculated values for the given percentiles
percentiles: list or numpy.ndarray
Percentiles to use for CDF matching

Returns
-------
uniq_perc_values: numpy.ndarray
Unique percentile values generated through fitting
the CDF of a beta distribution.

Raises
------
RuntimeError
If no fit could be found.
"""

# normalize between 0 and 1
min_value = np.min(perc_values)
perc_values = perc_values - min_value
max_value = np.max(perc_values)
perc_values = perc_values / max_value
percentiles = np.asanyarray(percentiles)
percentiles = percentiles / 100.0

p, ier = sc_opt.curve_fit(betainc,
percentiles,
perc_values)
uniq_perc_values = sc_special.betainc(p[0], p[1], percentiles)
uniq_perc_values = uniq_perc_values * max_value + min_value
return uniq_perc_values

[docs]def betainc(x, a, b):
return sc_special.betainc(a, b, x)

[docs]def element_iterable(el):
"""
Test if a element is iterable

Parameters
----------
el: object

Returns
-------
iterable: boolean
if True then then el is iterable
if Fales then not
"""
try:
el[0]
iterable = True
except (TypeError, IndexError):
iterable = False

return iterable

[docs]def ensure_iterable(el):
"""
Ensure that an object is iterable by putting it into a list.
Strings are handled slightly differently. They are
technically iterable but we want to keep the whole.

Parameters
----------
el: object

Returns
-------
iterable: list
[el]
"""
if type(el) == str:
return [el]
if not element_iterable(el):
return [el]
else:
return el