pyerrors.misc

  1import platform
  2import numpy as np
  3import scipy
  4import matplotlib
  5import matplotlib.pyplot as plt
  6import pandas as pd
  7import pickle
  8from .obs import Obs
  9from .version import __version__
 10
 11
 12def print_config():
 13    """Print information about version of python, pyerrors and dependencies."""
 14    config = {"system": platform.system(),
 15              "python": platform.python_version(),
 16              "pyerrors": __version__,
 17              "numpy": np.__version__,
 18              "scipy": scipy.__version__,
 19              "matplotlib": matplotlib.__version__,
 20              "pandas": pd.__version__}
 21
 22    for key, value in config.items():
 23        print(f"{key: <10}\t {value}")
 24
 25
 26def errorbar(x, y, axes=plt, **kwargs):
 27    """pyerrors wrapper for the errorbars method of matplotlib
 28
 29    Parameters
 30    ----------
 31    x : list
 32        A list of x-values which can be Obs.
 33    y : list
 34        A list of y-values which can be Obs.
 35    axes : (matplotlib.pyplot.axes)
 36        The axes to plot on. default is plt.
 37    """
 38    val = {}
 39    err = {}
 40    for name, comp in zip(["x", "y"], [x, y]):
 41        if all(isinstance(o, Obs) for o in comp):
 42            if not all(hasattr(o, 'e_dvalue') for o in comp):
 43                [o.gamma_method() for o in comp]
 44            val[name] = [o.value for o in comp]
 45            err[name] = [o.dvalue for o in comp]
 46        else:
 47            val[name] = comp
 48            err[name] = None
 49
 50        if f"{name}err" in kwargs:
 51            err[name] = kwargs.get(f"{name}err")
 52            kwargs.pop(f"{name}err", None)
 53
 54    axes.errorbar(val["x"], val["y"], xerr=err["x"], yerr=err["y"], **kwargs)
 55
 56
 57def dump_object(obj, name, **kwargs):
 58    """Dump object into pickle file.
 59
 60    Parameters
 61    ----------
 62    obj : object
 63        object to be saved in the pickle file
 64    name : str
 65        name of the file
 66    path : str
 67        specifies a custom path for the file (default '.')
 68
 69    Returns
 70    -------
 71    None
 72    """
 73    if 'path' in kwargs:
 74        file_name = kwargs.get('path') + '/' + name + '.p'
 75    else:
 76        file_name = name + '.p'
 77    with open(file_name, 'wb') as fb:
 78        pickle.dump(obj, fb)
 79
 80
 81def load_object(path):
 82    """Load object from pickle file.
 83
 84    Parameters
 85    ----------
 86    path : str
 87        path to the file
 88
 89    Returns
 90    -------
 91    object : Obs
 92        Loaded Object
 93    """
 94    with open(path, 'rb') as file:
 95        return pickle.load(file)
 96
 97
 98def pseudo_Obs(value, dvalue, name, samples=1000):
 99    """Generate an Obs object with given value, dvalue and name for test purposes
100
101    Parameters
102    ----------
103    value : float
104        central value of the Obs to be generated.
105    dvalue : float
106        error of the Obs to be generated.
107    name : str
108        name of the ensemble for which the Obs is to be generated.
109    samples: int
110        number of samples for the Obs (default 1000).
111
112    Returns
113    -------
114    res : Obs
115        Generated Observable
116    """
117    if dvalue <= 0.0:
118        return Obs([np.zeros(samples) + value], [name])
119    else:
120        for _ in range(100):
121            deltas = [np.random.normal(0.0, dvalue * np.sqrt(samples), samples)]
122            deltas -= np.mean(deltas)
123            deltas *= dvalue / np.sqrt((np.var(deltas) / samples)) / np.sqrt(1 + 3 / samples)
124            deltas += value
125            res = Obs(deltas, [name])
126            res.gamma_method(S=2, tau_exp=0)
127            if abs(res.dvalue - dvalue) < 1e-10 * dvalue:
128                break
129
130        res._value = float(value)
131
132        return res
133
134
135def gen_correlated_data(means, cov, name, tau=0.5, samples=1000):
136    """ Generate observables with given covariance and autocorrelation times.
137
138    Parameters
139    ----------
140    means : list
141        list containing the mean value of each observable.
142    cov : numpy.ndarray
143        covariance matrix for the data to be generated.
144    name : str
145        ensemble name for the data to be geneated.
146    tau : float or list
147        can either be a real number or a list with an entry for
148        every dataset.
149    samples : int
150        number of samples to be generated for each observable.
151
152    Returns
153    -------
154    corr_obs : list[Obs]
155        Generated observable list
156    """
157
158    assert len(means) == cov.shape[-1]
159    tau = np.asarray(tau)
160    if np.min(tau) < 0.5:
161        raise Exception('All integrated autocorrelations have to be >= 0.5.')
162
163    a = (2 * tau - 1) / (2 * tau + 1)
164    rand = np.random.multivariate_normal(np.zeros_like(means), cov * samples, samples)
165
166    # Normalize samples such that sample variance matches input
167    norm = np.array([np.var(o, ddof=1) / samples for o in rand.T])
168    rand = rand @ np.diag(np.sqrt(np.diag(cov))) @ np.diag(1 / np.sqrt(norm))
169
170    data = [rand[0]]
171    for i in range(1, samples):
172        data.append(np.sqrt(1 - a ** 2) * rand[i] + a * data[-1])
173    corr_data = np.array(data) - np.mean(data, axis=0) + means
174    return [Obs([dat], [name]) for dat in corr_data.T]
175
176
177def _assert_equal_properties(ol, otype=Obs):
178    otype = type(ol[0])
179    for o in ol[1:]:
180        if not isinstance(o, otype):
181            raise Exception("Wrong data type in list.")
182        for attr in ["reweighted", "e_content", "idl"]:
183            if hasattr(ol[0], attr):
184                if not getattr(ol[0], attr) == getattr(o, attr):
185                    raise Exception(f"All Obs in list have to have the same state '{attr}'.")