text:
"""
Distance statistics for planar point patterns

"""
__author__ = "Serge Rey sjsrey@gmail.com"
__all__ = ['DStatistic', 'G', 'F', 'J', 'K', 'L', 'Envelopes', 'Genv', 'Fenv', 'Jenv', 'Kenv', 'Lenv']

from .process import PoissonPointProcess as csr
import numpy as np
from matplotlib import pyplot as plt

class DStatistic(object):
    """
    Abstract Base Class for distance statistics.

    Parameters
    ----------
    name       : string
                 Name of the function. ("G", "F", "J", "K" or "L")

    Attributes
    ----------
    d          : array
                 The distance domain sequence.

    """
    def __init__(self, name):
        self.name = name

    def plot(self, qq=False):
        """
        Plot the distance function

        Parameters
        ----------
        qq: Boolean
            If False the statistic is plotted against distance. If Frue, the
            quantile-quantile plot is generated, observed vs. CSR.
        """

        # assuming mpl
        x = self.d
        if qq:
            plt.plot(self.ev, self._stat)
            plt.plot(self.ev, self.ev)
        else:
            plt.plot(x, self._stat, label='{}'.format(self.name))
            plt.ylabel("{}(d)".format(self.name))
            plt.xlabel('d')
            plt.plot(x, self.ev, label='CSR')
            plt.title("{} distance function".format(self.name))

class G(DStatistic):
    """
    Estimates the nearest neighbor distance distribution function G for a
    point pattern.

    Parameters
    ----------
    pp         : :class:`.PointPattern`
                 Point Pattern instance.
    intervals  : int
                 The length of distance domain sequence.
    dmin       : float
                 The minimum of the distance domain.
    dmax       : float
                 The maximum of the distance domain.
    d          : sequence
                 The distance domain sequence.
                 If d is specified, intervals, dmin and dmax are ignored.

    Attributes
    ----------
    name       : string
                 Name of the function. ("G", "F", "J", "K" or "L")
    d          : array
                 The distance domain sequence.
    G          : array
                 The cumulative nearest neighbor distance distribution over d.

    Notes
    -----
    In the analysis of planar point processes, the estimate of :math:`G` is
    typically compared to the value expected from a completely spatial
    random (CSR) process given as:

    .. math::

            G(d) = 1 - e^{-\lambda \pi  d^2}

    where :math:`\lambda` is the intensity (points per unit area) of the point
    process and :math:`d` is distance.

    For a clustered pattern, the empirical function will be above the
    expectation, while for a uniform pattern the empirical function falls below
    the expectation.

    """

    def __init__(self, pp, intervals=10, dmin=0.0, dmax=None, d=None):
        res = _g(pp, intervals, dmin, dmax, d)
        self.d = res[:, 0]
        self.G = self._stat = res[:, 1]
        self.ev = 1 - np.exp(-pp.lambda_window * np.pi * self.d * self.d)
        self.pp = pp
        super(G, self).__init__(name="G")

class F(DStatistic):
    """
    Estimates the empty space distribution function for a point pattern: F(d).

    Parameters
    ----------
    pp         : :class:`.PointPattern`
                 Point Pattern instance.
    n          : int
                 Number of empty space points (random points).
    intervals  : int
                 The length of distance domain sequence.
    dmin       : float
                 The minimum of the distance domain.
    dmax       : float
                 The maximum of the distance domain.
    d          : sequence
                 The distance domain sequence.
                 If d is specified, intervals, dmin and dmax are ignored.

    Attributes
    ----------
    d          : array
                 The distance domain sequence.
    G          : array
                 The cumulative empty space nearest event distance distribution
                 over d.

    Notes
    -----
    In the analysis of planar point processes, the estimate of :math:`F` is
    typically compared to the value expected from a process that displays
    complete spatial randomness (CSR):

    .. math::

            F(d) = 1 - e^{-\lambda \pi  d^2}

    where :math:`\lambda` is the intensity (points per unit area) of the point
    process and :math:`d` is distance.

    The expectation is identical to the expectation for the :class:`G` function
    for a CSR process.  However, for a clustered pattern, the empirical G
    function will be below the expectation, while for a uniform pattern the
    empirical function falls above the expectation.

    """

    def __init__(self, pp, n=100, intervals=10, dmin=0.0, dmax=None, d=None):
        res = _f(pp, n, intervals, dmin, dmax, d)
        self.d = res[:, 0]
        self.F = self._stat = res[:, 1]
        self.ev = 1 - np.exp(-pp.lambda_window * np.pi * self.d * self.d)
        super(F, self).__init__(name="F")

class J(DStatistic):
    """
    Estimates the J function for a point pattern :cite:`VanLieshout1996`

    Parameters
    ----------
    pp         : :class:`.PointPattern`
                 Point Pattern instance.
    n          : int
                 Number of empty space points (random points).
    intervals  : int
                 The length of distance domain sequence.
    dmin       : float
                 The minimum of the distance domain.
    dmax       : float
                 The maximum of the distance domain.
    d          : sequence
                 The distance domain sequence.
                 If d is specified, intervals, dmin and dmax are ignored.

    Attributes
    ----------
    d          : array
                 The distance domain sequence.
    j          : array
                 F function over d.

    Notes
    -----

    The :math:`J` function is a ratio of the hazard functions defined for
    :math:`G` and :math:`F`:

    .. math::

            J(d) = \\frac{1-G(d) }{1-F(d)}

    where :math:`G(d)` is the nearest neighbor distance distribution function
    (see :class:`G`)
    and :math:`F(d)` is the empty space function (see :class:`F`).

    For a CSR process the J function equals 1. Empirical values larger than 1
    are indicative of uniformity, while values below 1 suggest clustering.

    """
    def __init__(self, pp, n=100, intervals=10, dmin=0.0, dmax=None, d=None):
        res = _j(pp, n, intervals, dmin, dmax, d)
        self.d = res[:, 0]
        self.j = self._stat = res[:, 1]
        self.ev = self.j / self.j
        super(J, self).__init__(name="J")

class K(DStatistic):
    """
    Estimates the  K function for a point pattern.

    Parameters
    ----------
    pp         : :class:`.PointPattern`
                 Point Pattern instance.
    intervals  : int
                 The length of distance domain sequence.
    dmin       : float
                 The minimum of the distance domain.
    dmax       : float
                 The maximum of the distance domain.
    d          : sequence
                 The distance domain sequence.
                 If d is specified, intervals, dmin and dmax are ignored.

    Attributes
    ----------
    d          : array
                 The distance domain sequence.
    j          : array
                 K function over d.

    """
    def __init__(self, pp, intervals=10, dmin=0.0, dmax=None, d=None):
        res = _k(pp, intervals, dmin, dmax, d)
        self.d = res[:, 0]
        self.k = self._stat = res[:, 1]
        self.ev = np.pi * self.d * self.d
        super(K, self).__init__(name="K")

class L(DStatistic):
    """
    Estimates the l function for a point pattern.

    Parameters
    ----------
    pp         : :class:`.PointPattern`
                 Point Pattern instance.
    intervals  : int
                 The length of distance domain sequence.
    dmin       : float
                 The minimum of the distance domain.
    dmax       : float
                 The maximum of the distance domain.
    d          : sequence
                 The distance domain sequence.
                 If d is specified, intervals, dmin and dmax are ignored.

    Attributes
    ----------
    d          : array
                 The distance domain sequence.
    l          : array
                 L function over d.
    """
    def __init__(self, pp, intervals=10, dmin=0.0, dmax=None, d=None):
        res = _l(pp, intervals, dmin, dmax, d)
        self.d = res[:, 0]
        self.l = self._stat = res[:, 1]
        super(L, self).__init__(name="L")

    def plot(self):
        # assuming mpl
        x = self.d
        plt.plot(x, self._stat, label='{}'.format(self.name))
        plt.ylabel("{}(d)".format(self.name))
        plt.xlabel('d')
        plt.title("{} distance function".format(self.name))

def _g(pp, intervals=10, dmin=0.0, dmax=None, d=None):
    """
    Estimate the nearest neighbor distances function G.

    Parameters
    ----------
    pp       : :class:`.PointPattern`
               Point Pattern instance.
    intevals : int
               Number of intervals to evaluate F over.
    dmin     : float
               Lower limit of distance range.
    dmax     : float
               Upper limit of distance range. If dmax is None, dmax will be set
               to maximum nearest neighor distance.
    d        : sequence
               The distance domain sequence. If d is specified, intervals, dmin
               and dmax are ignored.

    Returns
    -------
             : array
               A 2-dimensional numpy array of 2 columns. The first column is
               the distance domain sequence for the point pattern. The second
               column is the cumulative nearest neighbor distance distribution.

    Notes
    -----
    See :class:`G`.

    """
    if d is None:
        w = pp.max_nnd/intervals
        if dmax:
            w = dmax/intervals
        d = [w*i for i in range(intervals + 2)]
    cdf = [0] * len(d)
    for i, d_i in enumerate(d):
        smaller = [nndi for nndi in pp.nnd if nndi <= d_i]
        cdf[i] = len(smaller)*1./pp.n
    return np.vstack((d, cdf)).T

def _f(pp, n=100, intervals=10, dmin=0.0, dmax=None, d=None):
    """
    F empty space function.

    Parameters
    ----------
    pp       : :class:`.PointPattern`
               Point Pattern instance.
    n        : int
               Number of empty space points (random points).
    intevals : int
               Number of intervals to evaluate F over.
    dmin     : float
               Lower limit of distance range.
    dmax     : float
               Upper limit of distance range. If dmax is None, dmax will be set
               to maximum nearest neighor distance.
    d        : sequence
               The distance domain sequence. If d is specified, intervals, dmin
               and dmax are ignored.

    Returns
    -------
             : array
               A 2-dimensional numpy array of 2 columns. The first column is
               the distance domain sequence for the point pattern. The second
               column is corresponding F function.

    Notes
    -----
    See :class:`.F`

    """

    # get a csr pattern in window of pp
    c = csr(pp.window, n, 1, asPP=True).realizations[0]
    # for each point in csr pattern find the closest point in pp and the
    # associated distance
    nnids, nnds = pp.knn_other(c, k=1)

    if d is None:
        w = pp.max_nnd/intervals
        if dmax:
            w = dmax/intervals
        d = [w*i for i in range(intervals + 2)]
    cdf = [0] * len(d)

    for i, d_i in enumerate(d):
        smaller = [nndi for nndi in nnds if nndi <= d_i]
        cdf[i] = len(smaller)*1./n
    return np.vstack((d, cdf)).T

def _j(pp, n=100, intervals=10, dmin=0.0, dmax=None, d=None):
    """
    J function: Ratio of hazard functions for F and G.

    Parameters
    ----------
    pp       : :class:`.PointPattern`
               Point Pattern instance.
    n        : int
               Number of empty space points (random points).
    intevals : int
               Number of intervals to evaluate F over.
    dmin     : float
               Lower limit of distance range.
    dmax     : float
               Upper limit of distance range. If dmax is None, dmax will be set
               to maximum nearest neighor distance.
    d        : sequence
               The distance domain sequence. If d is specified, intervals, dmin
               and dmax are ignored.

    Returns
    -------
             : array
               A 2-dimensional numpy array of 2 columns. The first column is
               the distance domain sequence for the point pattern. The second
               column is corresponding J function.

    Notes
    -----
    See :class:`.J`

    """

    F = _f(pp, n, intervals=intervals, dmin=dmin, dmax=dmax, d=d)
    G = _g(pp, intervals=intervals, dmin=dmin, dmax=dmax, d=d)
    FC = 1 - F[:, 1]
    GC = 1 - G[:, 1]
    last_id = len(GC) + 1
    if np.any(FC == 0):
        last_id = np.where(FC == 0)[0][0]

    return np.vstack((F[:last_id, 0], GC[:last_id]/FC[:last_id])).T

def _k(pp, intervals=10, dmin=0.0, dmax=None, d=None):
    """
    Interevent K function.

    Parameters
    ----------
    pp       : :class:`.PointPattern`
               Point Pattern instance.
    n        : int
               Number of empty space points (random points).
    intevals : int
               Number of intervals to evaluate F over.
    dmin     : float
               Lower limit of distance range.
    dmax     : float
               Upper limit of distance range. If dmax is None, dmax will be set
               to length of bounding box diagonal.
    d        : sequence
               The distance domain sequence. If d is specified, intervals, dmin
               and dmax are ignored.

    Returns
    -------
    kcdf     : array
               A 2-dimensional numpy array of 2 columns. The first column is
               the distance domain sequence for the point pattern. The second
               column is corresponding K function.

    Notes
    -----
    See :class:`.K`

    """

    if d is None:
        # use length of bounding box diagonal as max distance
        bb = pp.mbb
        dbb = np.sqrt((bb[0]-bb[2])**2 + (bb[1]-bb[3])**2)
        w = dbb/intervals
        if dmax:
            w = dmax/intervals
        d = [w*i for i in range(intervals + 2)]
    den = pp.lambda_window * pp.n * 2.
    kcdf = np.asarray([(di, len(pp.tree.query_pairs(di))/den) for di in d])
    return kcdf

def _l(pp, intervals=10, dmin=0.0, dmax=None, d=None):
    """
    Interevent L function.

    Parameters
    ----------
    pp       : :class:`.PointPattern`
               Point Pattern instance.
    n        : int
               Number of empty space points (random points).
    intevals : int
               Number of intervals to evaluate F over.
    dmin     : float
               Lower limit of distance range.
    dmax     : float
               Upper limit of distance range. If dmax is None, dmax will be set
               to length of bounding box diagonal.
    d        : sequence
               The distance domain sequence. If d is specified, intervals, dmin
               and dmax are ignored.

    Returns
    -------
    kf       : array
               A 2-dimensional numpy array of 2 columns. The first column is
               the distance domain sequence for the point pattern. The second
               column is corresponding L function.

    Notes
    -----
    See :class:`.L`

    """

    kf = _k(pp, intervals, dmin, dmax, d)
    kf[:, 1] = np.sqrt(kf[:, 1] / np.pi) - kf[:, 0]
    return kf

class Envelopes(object):
    """
    Abstract base class for simulation envelopes.

    Parameters
    ----------
    pp          : :class:`.PointPattern`
                  Point Pattern instance.
    intervals   : int
                  The length of distance domain sequence. Default is 10.
    dmin        : float
                  The minimum of the distance domain.
    dmax        : float
                  The maximum of the distance domain.
    d           : sequence
                  The distance domain sequence.
                  If d is specified, intervals, dmin and dmax are ignored.
    pct         : float
                  1-alpha, alpha is the significance level. Default is 0.05,
                  1-alpha is the confidence level for the envelope.
    realizations: :class:`.PointProcess`
                  Point process instance with more than 1 realizations.

    Attributes
    ----------
    name        : string
                  Name of the function. ("G", "F", "J", "K" or "L")
    observed    : array
                  A 2-dimensional numpy array of 2 columns. The first column is
                  the distance domain sequence for the observed point pattern.
                  The second column is the specific function ("G", "F", "J",
                  "K" or "L") over the distance domain sequence for the
                  observed point pattern.
    low         : array
                  A 1-dimensional numpy array. Lower bound of the simulation
                  envelope.
    high        : array
                  A 1-dimensional numpy array. Higher bound of the simulation
                  envelope.
    mean        : array
                  A 1-dimensional numpy array. Mean values of the simulation
                  envelope.

    """
    def __init__(self, *args,  **kwargs):
        # setup arguments
        self.name = kwargs['name']

        # calculate observed function
        self.pp = args[0]
        self.observed = self.calc(*args, **kwargs)
        self.d = self.observed[:, 0]  # domain to be used in all realizations

        # do realizations
        self.mapper(kwargs['realizations'])

    def mapper(self, realizations):
        reals = realizations.realizations
        res = np.asarray([self.calc(reals[p]) for p in reals])

        # When calculating the J function for all the simulations, the length
        # of the returned interval domains might be different.

        if self.name == "J":
            res = []
            for p in reals:
                j = self.calc(reals[p])
                if j.shape[0] < self.d.shape[0]:
                    diff = self.d.shape[0]-j.shape[0]
                    for i in range(diff):
                        j = np.append(j, [[self.d[i+diff], np.inf]], axis=0)
                res.append(j)
            res = np.array(res)

        res = res[:, :, -1]
        res.sort(axis=0)
        nres = len(res)
        self.low = res[np.int(nres * self.pct/2.)]
        self.high = res[np.int(nres * (1-self.pct/2.))]
        self.mean = res.mean(axis=0)

    def calc(self, *args, **kwargs):
        print('implement in subclass')

    def plot(self):
        # assuming mpl
        x = self.d
        plt.plot(x, self.observed[:, 1], label='{}'.format(self.name))
        plt.plot(x, self.mean, 'g-.', label='CSR')
        plt.plot(x, self.low, 'r-.', label='LB')
        plt.plot(x, self.high, 'r-.', label="UB")
        plt.ylabel("{}(d)".format(self.name))
        plt.xlabel('d')
        plt.title("{} Simulation Envelopes".format(self.name))
        plt.legend(loc=0)

class Genv(Envelopes):
    """
    Simulation envelope for G function.

    Parameters
    ----------
    pp          : :class:`.PointPattern`
                  Point Pattern instance.
    intervals   : int
                  The length of distance domain sequence. Default is 10.
    dmin        : float
                  The minimum of the distance domain.
    dmax        : float
                  Upper limit of distance range. If dmax is None, dmax will be
                  set to maximum nearest neighbor distance.
    d           : sequence
                  The distance domain sequence.
                  If d is specified, intervals, dmin and dmax are ignored.
    pct         : float
                  1-alpha, alpha is the significance level. Default is 0.05,
                  which means 95% confidence level for the envelopes.
    realizations: :class:`.PointProcess`
                  Point process instance with more than 1 realizations.

    Attributes
    ----------
    name        : string
                  Name of the function. ("G", "F", "J", "K" or "L")
    observed    : array
                  A 2-dimensional numpy array of 2 columns. The first column is
                  the distance domain sequence for the observed point pattern.
                  The second column is cumulative nearest neighbor distance
                  distribution (G function) for the observed point pattern.
    low         : array
                  A 1-dimensional numpy array. Lower bound of the simulation
                  envelope.
    high        : array
                  A 1-dimensional numpy array. Higher bound of the simulation
                  envelope.
    mean        : array
                  A 1-dimensional numpy array. Mean values of the simulation
                  envelope.

    Examples
    --------
    .. plot::

       >>> import pysal.lib as ps
       >>> from pointpats import Genv, PoissonPointProcess, Window
       >>> from pysal.lib.cg import shapely_ext
       >>> va = ps.io.open(ps.examples.get_path("vautm17n.shp"))
       >>> polys = [shp for shp in va]
       >>> state = shapely_ext.cascaded_union(polys)
       >>> pp = PoissonPointProcess(Window(state.parts), 100, 1, asPP=True).realizations[0]
       >>> csrs = PoissonPointProcess(pp.window, 100, 100, asPP=True)
       >>> genv_bb = Genv(pp, realizations=csrs)
       >>> genv_bb.plot()

    """

    def __init__(self, pp, intervals=10, dmin=0.0, dmax=None, d=None, pct=0.05,
                 realizations=None):
        self.pp = pp
        self.intervals = intervals
        self.dmin = dmin
        self.dmax = dmax
        self.d = d
        self.pct = pct
        super(Genv, self).__init__(pp, realizations=realizations, name="G")

    def calc(self, *args, **kwargs):
        pp = args[0]
        return _g(pp, intervals=self.intervals, dmin=self.dmin, dmax=self.dmax,
                  d=self.d)

class Fenv(Envelopes):
    """
    Simulation envelope for F function.

    Parameters
    ----------
    pp          : :class:`.PointPattern`
                  Point Pattern instance.
    n           : int
                  Number of empty space points (random points).
    intervals   : int
                  The length of distance domain sequence. Default is 10.
    dmin        : float
                  The minimum of the distance domain.
    dmax        : float
                  Upper limit of distance range. If dmax is None, dmax will be
                  set to maximum nearest neighbor distance.
    d           : sequence
                  The distance domain sequence.
                  If d is specified, intervals, dmin and dmax are ignored.
    pct         : float
                  1-alpha, alpha is the significance level. Default is 0.05,
                  which means 95% confidence level for the envelopes.
    realizations: :class:`.PointProcess`
                  Point process instance with more than 1 realizations.

    Attributes
    ----------
    name        : string
                  Name of the function. ("G", "F", "J", "K" or "L")
    observed    : array
                  A 2-dimensional numpy array of 2 columns. The first column is
                  the distance domain sequence for the observed point pattern.
                  The second column is F function for the observed point
                  pattern.
    low         : array
                  A 1-dimensional numpy array. Lower bound of the simulation
                  envelope.
    high        : array
                  A 1-dimensional numpy array. Higher bound of the simulation
                  envelope.
    mean        : array
                  A 1-dimensional numpy array. Mean values of the simulation
                  envelope.

    Examples
    --------
    .. plot::

       >>> import pysal.lib as ps
       >>> from pysal.lib.cg import shapely_ext
       >>> from pointpats import PoissonPointProcess,Window,Fenv
       >>> va = ps.io.open(ps.examples.get_path("vautm17n.shp"))
       >>> polys = [shp for shp in va]
       >>> state = shapely_ext.cascaded_union(polys)
       >>> pp = PoissonPointProcess(Window(state.parts), 100, 1, asPP=True).realizations[0]
       >>> csrs = PoissonPointProcess(pp.window, 100, 100, asPP=True)
       >>> fenv = Fenv(pp, realizations=csrs)
       >>> fenv.plot()

    """
    def __init__(self, pp, n=100, intervals=10, dmin=0.0, dmax=None, d=None,
                 pct=0.05, realizations=None):
        self.pp = pp
        self.n = n
        self.intervals = intervals
        self.dmin = dmin
        self.dmax = dmax
        self.d = d
        self.pct = pct
        super(Fenv, self).__init__(pp, realizations=realizations, name="F")

    def calc(self, *args, **kwargs):
        pp = args[0]
        return _f(pp, self.n, intervals=self.intervals, dmin=self.dmin,
                  dmax=self.dmax, d=self.d)

class Jenv(Envelopes):
    """
    Simulation envelope for J function.

    Parameters
    ----------
    pp          : :class:`.PointPattern`
                  Point Pattern instance.
    n           : int
                  Number of empty space points (random points).
    intervals   : int
                  The length of distance domain sequence. Default is 10.
    dmin        : float
                  The minimum of the distance domain.
    dmax        : float
                  Upper limit of distance range. If dmax is None, dmax will be
                  set to maximum nearest neighbor distance.
    d           : sequence
                  The distance domain sequence.
                  If d is specified, intervals, dmin and dmax are ignored.
    pct         : float
                  1-alpha, alpha is the significance level. Default is 0.05,
                  which means 95% confidence level for the envelopes.
    realizations: :class:`.PointProcess`
                  Point process instance with more than 1 realizations.

    Attributes
    ----------
    name        : string
                  Name of the function. ("G", "F", "J", "K" or "L")
    observed    : array
                  A 2-dimensional numpy array of 2 columns. The first column is
                  the distance domain sequence for the observed point pattern.
                  The second column is J function for the observed point
                  pattern.
    low         : array
                  A 1-dimensional numpy array. Lower bound of the simulation
                  envelope.
    high        : array
                  A 1-dimensional numpy array. Higher bound of the simulation
                  envelope.
    mean        : array
                  A 1-dimensional numpy array. Mean values of the simulation
                  envelope.

    Examples
    --------
    .. plot::

       >>> import pysal.lib as ps
       >>> from pointpats import Jenv, PoissonPointProcess, Window
       >>> from pysal.lib.cg import shapely_ext
       >>> va = ps.io.open(ps.examples.get_path("vautm17n.shp"))
       >>> polys = [shp for shp in va]
       >>> state = shapely_ext.cascaded_union(polys)
       >>> pp = PoissonPointProcess(Window(state.parts), 100, 1, asPP=True).realizations[0]
       >>> csrs = PoissonPointProcess(pp.window, 100, 100, asPP=True)
       >>> jenv = Jenv(pp, realizations=csrs)
       >>> jenv.plot()

    """
    def __init__(self, pp, n=100, intervals=10, dmin=0.0, dmax=None, d=None,
                 pct=0.05, realizations=None):
        self.pp = pp
        self.n = n
        self.intervals = intervals
        self.dmin = dmin
        self.dmax = dmax
        self.d = d
        self.pct = pct
        super(Jenv, self).__init__(pp, realizations=realizations, name="J")

    def calc(self, *args, **kwargs):
        pp = args[0]
        return _j(pp, self.n, intervals=self.intervals, dmin=self.dmin,
                  dmax=self.dmax, d=self.d)

class Kenv(Envelopes):
    """
    Simulation envelope for K function.

    Parameters
    ----------
    pp          : :class:`.PointPattern`
                  Point Pattern instance.
    intervals   : int
                  The length of distance domain sequence. Default is 10.
    dmin        : float
                  The minimum of the distance domain.
    dmax        : float
                  Upper limit of distance range. If dmax is None, dmax will be
                  set to maximum nearest neighbor distance.
    d           : sequence
                  The distance domain sequence.
                  If d is specified, intervals, dmin and dmax are ignored.
    pct         : float
                  1-alpha, alpha is the significance level. Default is 0.05,
                  which means 95% confidence level for the envelope.
    realizations: :class:`.PointProcess`
                  Point process instance with more than 1 realizations.

    Attributes
    ----------
    name        : string
                  Name of the function. ("G", "F", "J", "K" or "L")
    observed    : array
                  A 2-dimensional numpy array of 2 columns. The first column is
                  the distance domain sequence for the observed point pattern.
                  The second column is K function for the observed point
                  pattern.
    low         : array
                  A 1-dimensional numpy array. Lower bound of the simulation
                  envelope.
    high        : array
                  A 1-dimensional numpy array. Higher bound of the simulation
                  envelope.
    mean        : array
                  A 1-dimensional numpy array. Mean values of the simulation
                  envelope.

    Examples
    --------
    .. plot::

       >>> import pysal.lib as ps
       >>> from pointpats import Kenv, PoissonPointProcess, Window
       >>> from pysal.lib.cg import shapely_ext
       >>> va = ps.io.open(ps.examples.get_path("vautm17n.shp"))
       >>> polys = [shp for shp in va]
       >>> state = shapely_ext.cascaded_union(polys)
       >>> pp = PoissonPointProcess(Window(state.parts), 100, 1, asPP=True).realizations[0]
       >>> csrs = PoissonPointProcess(pp.window, 100, 100, asPP=True)
       >>> kenv = Kenv(pp, realizations=csrs)
       >>> kenv.plot()

    """
    def __init__(self, pp, intervals=10, dmin=0.0, dmax=None, d=None,
                 pct=0.05, realizations=None):
        self.pp = pp
        self.intervals = intervals
        self.dmin = dmin
        self.dmax = dmax
        self.d = d
        self.pct = pct
        super(Kenv, self).__init__(pp, realizations=realizations, name="K")

    def calc(self, *args, **kwargs):
        pp = args[0]
        return _k(pp, intervals=self.intervals, dmin=self.dmin, dmax=self.dmax,
                  d=self.d)

class Lenv(Envelopes):
    """
    Simulation envelope for L function.

    Parameters
    ----------
    pp          : :class:`.PointPattern`
                  Point Pattern instance.
    intervals   : int
                  The length of distance domain sequence. Default is 10.
    dmin        : float
                  The minimum of the distance domain.
    dmax        : float
                  Upper limit of distance range. If dmax is None, dmax will be
                  set to maximum nearest neighbor distance.
    d           : sequence
                  The distance domain sequence.
                  If d is specified, intervals, dmin and dmax are ignored.
    pct         : float
                  1-alpha, alpha is the significance level. Default is 0.05,
                  which means 95% confidence level for the envelopes.
    realizations: :class:`.PointProcess`
                  Point process instance with more than 1 realizations.

    Attributes
    ----------
    name        : string
                  Name of the function. ("G", "F", "J", "K" or "L")
    observed    : array
                  A 2-dimensional numpy array of 2 columns. The first column is
                  the distance domain sequence for the observed point pattern.
                  The second column is L function for the observed point
                  pattern.
    low         : array
                  A 1-dimensional numpy array. Lower bound of the simulation
                  envelope.
    high        : array
                  A 1-dimensional numpy array. Higher bound of the simulation
                  envelope.
    mean        : array
                  A 1-dimensional numpy array. Mean values of the simulation
                  envelope.

    Examples
    --------
    .. plot::

       >>> import pysal.lib as ps
       >>> from pointpats import Lenv, PoissonPointProcess, Window
       >>> from pysal.lib.cg import shapely_ext
       >>> va = ps.io.open(ps.examples.get_path("vautm17n.shp"))
       >>> polys = [shp for shp in va]
       >>> state = shapely_ext.cascaded_union(polys)
       >>> pp = PoissonPointProcess(Window(state.parts), 100, 1, asPP=True).realizations[0]
       >>> csrs = PoissonPointProcess(pp.window, 100, 100, asPP=True)
       >>> lenv = Lenv(pp, realizations=csrs)
       >>> lenv.plot()

    """

    def __init__(self, pp, intervals=10, dmin=0.0, dmax=None, d=None,
                 pct=0.05, realizations=None):
        self.pp = pp
        self.intervals = intervals
        self.dmin = dmin
        self.dmax = dmax
        self.d = d
        self.pct = pct
        super(Lenv, self).__init__(pp, realizations=realizations, name="L")

    def calc(self, *args, **kwargs):
        pp = args[0]
        return _l(pp, intervals=self.intervals, dmin=self.dmin, dmax=self.dmax,
                  d=self.d)

max_stars_repo_name: ocefpaf/pysal

max_issues_repo_name: ocefpaf/pysal

max_forks_repo_name: ocefpaf/pysal


text:
#!/usr/bin/env python3
# Copyright (c) Facebook, Inc. and its affiliates.
#
# This source code is licensed under the MIT license found in the
# LICENSE file in the root directory of this source tree.

import numpy as np
import torch
from ax.exceptions.model import ModelError
from ax.models.torch.utils import (
    _generate_sobol_points,
    is_noiseless,
    normalize_indices,
    subset_model,
    tensor_callable_to_array_callable,
)
from ax.utils.common.testutils import TestCase
from botorch.models import HeteroskedasticSingleTaskGP, ModelListGP, SingleTaskGP
from torch import Tensor

class TorchUtilsTest(TestCase):
    def test_is_noiseless(self):
        x = torch.zeros(1, 1)
        y = torch.zeros(1, 1)
        se = torch.zeros(1, 1)
        model = SingleTaskGP(x, y)
        self.assertTrue(is_noiseless(model))
        model = HeteroskedasticSingleTaskGP(x, y, se)
        self.assertFalse(is_noiseless(model))
        with self.assertRaises(ModelError):
            is_noiseless(ModelListGP())

    def testNormalizeIndices(self):
        indices = [0, 2]
        nlzd_indices = normalize_indices(indices, 3)
        self.assertEqual(nlzd_indices, indices)
        nlzd_indices = normalize_indices(indices, 4)
        self.assertEqual(nlzd_indices, indices)
        indices = [0, -1]
        nlzd_indices = normalize_indices(indices, 3)
        self.assertEqual(nlzd_indices, [0, 2])
        with self.assertRaises(ValueError):
            nlzd_indices = normalize_indices([3], 3)
        with self.assertRaises(ValueError):
            nlzd_indices = normalize_indices([-4], 3)

    def testSubsetModel(self):
        x = torch.zeros(1, 1)
        y = torch.rand(1, 2)
        obj_t = torch.rand(2)
        model = SingleTaskGP(x, y)
        self.assertEqual(model.num_outputs, 2)
        # basic test, can subset
        obj_weights = torch.tensor([1.0, 0.0])
        subset_model_results = subset_model(model, obj_weights)
        model_sub = subset_model_results.model
        obj_weights_sub = subset_model_results.objective_weights
        ocs_sub = subset_model_results.outcome_constraints
        obj_t_sub = subset_model_results.objective_thresholds
        self.assertIsNone(ocs_sub)
        self.assertIsNone(obj_t_sub)
        self.assertEqual(model_sub.num_outputs, 1)
        self.assertTrue(torch.equal(obj_weights_sub, torch.tensor([1.0])))
        # basic test, cannot subset
        obj_weights = torch.tensor([1.0, 2.0])
        subset_model_results = subset_model(model, obj_weights)
        model_sub = subset_model_results.model
        obj_weights_sub = subset_model_results.objective_weights
        ocs_sub = subset_model_results.outcome_constraints
        obj_t_sub = subset_model_results.objective_thresholds
        self.assertIsNone(ocs_sub)
        self.assertIsNone(obj_t_sub)
        self.assertIs(model_sub, model)  # check identity
        self.assertIs(obj_weights_sub, obj_weights)  # check identity
        # test w/ outcome constraints, can subset
        obj_weights = torch.tensor([1.0, 0.0])
        ocs = (torch.tensor([[1.0, 0.0]]), torch.tensor([1.0]))
        subset_model_results = subset_model(model, obj_weights, ocs)
        model_sub = subset_model_results.model
        obj_weights_sub = subset_model_results.objective_weights
        ocs_sub = subset_model_results.outcome_constraints
        obj_t_sub = subset_model_results.objective_thresholds
        self.assertEqual(model_sub.num_outputs, 1)
        self.assertIsNone(obj_t_sub)
        self.assertTrue(torch.equal(obj_weights_sub, torch.tensor([1.0])))
        self.assertTrue(torch.equal(ocs_sub[0], torch.tensor([[1.0]])))
        self.assertTrue(torch.equal(ocs_sub[1], torch.tensor([1.0])))
        # test w/ outcome constraints, cannot subset
        obj_weights = torch.tensor([1.0, 0.0])
        ocs = (torch.tensor([[0.0, 1.0]]), torch.tensor([1.0]))
        subset_model_results = subset_model(model, obj_weights, ocs)
        model_sub = subset_model_results.model
        obj_weights_sub = subset_model_results.objective_weights
        ocs_sub = subset_model_results.outcome_constraints
        obj_t_sub = subset_model_results.objective_thresholds
        self.assertIs(model_sub, model)  # check identity
        self.assertIsNone(obj_t_sub)
        self.assertIs(obj_weights_sub, obj_weights)  # check identity
        self.assertIs(ocs_sub, ocs)  # check identity
        # test w/ objective thresholds, cannot subset
        obj_weights = torch.tensor([1.0, 0.0])
        ocs = (torch.tensor([[0.0, 1.0]]), torch.tensor([1.0]))
        subset_model_results = subset_model(model, obj_weights, ocs, obj_t)
        model_sub = subset_model_results.model
        obj_weights_sub = subset_model_results.objective_weights
        ocs_sub = subset_model_results.outcome_constraints
        obj_t_sub = subset_model_results.objective_thresholds
        self.assertIs(model_sub, model)  # check identity
        self.assertIs(obj_t, obj_t_sub)
        self.assertIs(obj_weights_sub, obj_weights)  # check identity
        self.assertTrue(torch.equal(subset_model_results.indices, torch.tensor([0, 1])))
        self.assertIs(ocs_sub, ocs)  # check identity
        # test w/ objective thresholds, can subset
        obj_weights = torch.tensor([1.0, 0.0])
        ocs = (torch.tensor([[1.0, 0.0]]), torch.tensor([1.0]))
        subset_model_results = subset_model(model, obj_weights, ocs, obj_t)
        model_sub = subset_model_results.model
        obj_weights_sub = subset_model_results.objective_weights
        ocs_sub = subset_model_results.outcome_constraints
        obj_t_sub = subset_model_results.objective_thresholds
        self.assertTrue(torch.equal(subset_model_results.indices, torch.tensor([0])))
        self.assertEqual(model_sub.num_outputs, 1)
        self.assertTrue(torch.equal(obj_weights_sub, torch.tensor([1.0])))
        self.assertTrue(torch.equal(obj_t_sub, obj_t[:1]))
        self.assertTrue(torch.equal(ocs_sub[0], torch.tensor([[1.0]])))
        self.assertTrue(torch.equal(ocs_sub[1], torch.tensor([1.0])))
        # test unsupported
        yvar = torch.ones(1, 2)
        model = HeteroskedasticSingleTaskGP(x, y, yvar)
        subset_model_results = subset_model(model, obj_weights)
        model_sub = subset_model_results.model
        obj_weights_sub = subset_model_results.objective_weights
        ocs_sub = subset_model_results.outcome_constraints
        obj_t_sub = subset_model_results.objective_thresholds
        self.assertIsNone(ocs_sub)
        self.assertIs(model_sub, model)  # check identity
        self.assertIs(obj_weights_sub, obj_weights)  # check identity
        # test error on size inconsistency
        obj_weights = torch.ones(3)
        with self.assertRaises(RuntimeError):
            subset_model(model, obj_weights)

    def testGenerateSobolPoints(self):
        bounds = [(0.0, 1.0) for _ in range(3)]
        linear_constraints = (
            torch.tensor([[1, -1, 0]], dtype=torch.double),
            torch.tensor([[0]], dtype=torch.double),
        )

        def test_rounding_func(x: Tensor) -> Tensor:
            return x

        gen_sobol = _generate_sobol_points(
            n_sobol=100,
            bounds=bounds,
            device=torch.device("cpu"),
            linear_constraints=linear_constraints,
            rounding_func=test_rounding_func,
        )
        self.assertEqual(len(gen_sobol), 100)
        self.assertIsInstance(gen_sobol, Tensor)

    def testTensorCallableToArrayCallable(self):
        def tensor_func(x: Tensor) -> Tensor:
            return np.exp(x)

        new_func = tensor_callable_to_array_callable(
            tensor_func=tensor_func, device=torch.device("cpu")
        )
        self.assertTrue(callable(new_func))
        self.assertIsInstance(new_func(np.array([1.0, 2.0])), np.ndarray)

ext: py

max_stars_repo_name: saitcakmak/Ax

max_issues_repo_name: saitcakmak/Ax

max_forks_repo_name: saitcakmak/Ax

reason: import numpy


text:
from ...isa.inst import *
import numpy as np
import math

class Vssex_v(Inst):
    name = 'vssex.v'

    def golden(self):

        if 'isExcept' in self:
            if self['isExcept'] > 0:
                return np.zeros(self['vl'], dtype=self['vs3'].dtype)
        vl = self['vl']
        emul = self['eew'] /self['sew']*factor_lmul[self['lmul']]
        emul = 1 if emul < 1 else int(emul)
        vlmax = int(emul * self['VLEN'] // (self['vs3'].itemsize*8))
        eew = self['eew']//8
        stride = self['rs2']//eew

        vs3 = self['vs3']

        if 'start' in self:
            start = self['start']
        else:
            start = 0

        if 'origin' in self:
            origin = self['origin']
        else:
            origin = np.zeros(vl, dtype=self['vs3'].dtype)

        if 'mask' in self:
            mask = np.unpackbits(self['mask'], bitorder='little')[0: self['vl']]
        else :
            mask = np.ones(self['vl'], dtype=np.uint8)

        if stride == 0:
            tmp = 0
            for i in range(start, vl):
                if mask[i] != 0:
                    tmp = vs3[i]
            origin[0: vl] = tmp
        else:
            origin[start: vl] = self.masked(vs3[0:vl], origin[0: vl])[start: vl]

        return origin

max_stars_repo_name: ultrafive/riscv-pvp

max_issues_repo_name: ultrafive/riscv-pvp

max_forks_repo_name: ultrafive/riscv-pvp


text:
from CHECLabPy.plotting.setup import Plotter
import numpy as np
import matplotlib as mpl
from matplotlib.colors import LogNorm
from matplotlib.collections import PatchCollection
from matplotlib.patches import Rectangle, Circle
from CHECLabPy.utils.mapping import get_clp_mapping_from_tc_mapping
import os
from copy import copy

class CameraImage(Plotter):
    def __init__(self, xpix, ypix, size, cmap=None, **kwargs):
        """
        Create a camera-image plot

        Parameters
        ----------
        xpix : ndarray
            The X positions of the pixels/superpixels/TMs
        ypix : ndarray
            The Y positions of the pixels/superpixels/TMs
        size : float
            The size of the pixels/superpixels/TMs
        kwargs
            Arguments passed to `CHECLabPy.plottong.setup.Plotter`
        """
        super().__init__(**kwargs)

        rc = {
            "xtick.direction": 'out',
            "ytick.direction": 'out',
        }
        mpl.rcParams.update(rc)

        self._image = None
        self._mapping = None
        self.colorbar = None
        self.autoscale = True

        self.xpix = xpix
        self.ypix = ypix

        assert self.xpix.size == self.ypix.size
        self.n_pixels = self.xpix.size

        patches = []
        for xx, yy in zip(self.xpix, self.ypix):
            rr = size + 0.0001  # extra size to pixels to avoid aliasing
            poly = Rectangle(
                (xx - rr / 2., yy - rr / 2.),
                width=rr,
                height=rr,
                fill=True,
            )
            patches.append(poly)

        self.pixels = PatchCollection(patches, linewidth=0, cmap=cmap)
        self.ax.add_collection(self.pixels)
        self.pixels.set_array(np.zeros(self.n_pixels))

        self.ax.set_aspect('equal', 'datalim')
        self.ax.set_xlabel("X position (m)")
        self.ax.set_ylabel("Y position (m)")
        self.ax.autoscale_view()
        self.ax.axis('off')

    @property
    def image(self):
        return self._image

    @image.setter
    def image(self, val):
        assert val.size == self.n_pixels

        self._image = val

        self.pixels.set_array(np.ma.masked_invalid(val))
        self.pixels.changed()
        if self.autoscale:
            self.pixels.autoscale()  # Updates the colorbar
        self.ax.figure.canvas.draw()

    def save(self, output_path, **kwargs):
        super().save(output_path, **kwargs)
        if output_path.endswith('.pdf'):
            try:
                self.crop(output_path)
            except ModuleNotFoundError:
                pass

    @staticmethod
    def crop(path):
        from PyPDF2 import PdfFileWriter, PdfFileReader
        with open(path, "rb") as in_f:
            input1 = PdfFileReader(in_f)
            output = PdfFileWriter()

            num_pages = input1.getNumPages()
            for i in range(num_pages):
                page = input1.getPage(i)
                page.cropBox.lowerLeft = (100, 20)
                page.cropBox.upperRight = (340, 220)
                output.addPage(page)

            pdf_path = os.path.splitext(path)[0] + "_cropped.pdf"
            with open(pdf_path, "wb") as out_f:
                output.write(out_f)
                print("Cropped figure saved to: {}".format(pdf_path))

    def set_cmap(self,cmap="viridis"):
        self.pixels.set_cmap(cmap)

    def add_colorbar(self, label='', pad=-0.2, ax=None,**kwargs):
        if ax is None:
            ax = self.ax
        self.colorbar = self.ax.figure.colorbar(
            self.pixels, label=label, pad=pad, ax=ax,**kwargs
        )

    def set_limits_minmax(self, zmin, zmax):
        """
        Set the color scale limits from min to max
        """
        self.pixels.set_clim(zmin, zmax)
        self.autoscale = False

    def set_z_log(self):
        self.pixels.norm = LogNorm()
        self.pixels.autoscale()

    def reset_limits(self):
        """
        Reset to auto color scale limits
        """
        self.autoscale = True
        self.pixels.autoscale()

    def annotate_on_telescope_up(self):
        """
        Add an arrow indicating where "ON-Telescope-UP" is
        """
        if self._mapping is not None:
            axl = self._mapping.metadata['fOTUpX_l']
            ayl = self._mapping.metadata['fOTUpY_l']
            adx = self._mapping.metadata['fOTUpX_u'] - axl
            ady = self._mapping.metadata['fOTUpY_u'] - ayl
            text = "ON-Telescope UP"
            self.ax.arrow(axl, ayl, adx, ady, head_width=0.01,
                          head_length=0.01, fc='r', ec='r')
            self.ax.text(axl, ayl, text, fontsize=4, color='r',
                         ha='center', va='bottom')
        else:
            print("Cannot annotate, no mapping attached to class")

    def add_text_to_pixel(self, pixel, value, size=3, color='w', **kwargs):
        """
        Add a text label to a single pixel

        Parameters
        ----------
        pixel : int
        value : str or float
        size : int
            Font size
        color : str
            Color of the text
        kwargs
            Named arguments to pass to matplotlib.pyplot.text
        """
        pos_x = self.xpix[pixel]
        pos_y = self.ypix[pixel]
        self.ax.text(pos_x, pos_y, value,
                     fontsize=size, color=color, ha='center', **kwargs)

    def add_pixel_text(self, values, size=3, color='w', **kwargs):
        """
        Add a text label to each pixel

        Parameters
        ----------
        values : ndarray
        size : int
            Font size
        color : str
            Color of the text
        kwargs
            Named arguments to pass to matplotlib.pyplot.text
        """
        assert values.size == self.n_pixels
        for pixel in range(self.n_pixels):
            self.add_text_to_pixel(pixel, values[pixel], size, color, **kwargs)

    def highlight_pixels(self, pixels, color='g', linewidth=0.5, alpha=0.75):
        """
        Highlight the given pixels with a colored line around them

        Parameters
        ----------
        pixels : index-like
            The pixels to highlight.
            Can either be a list or array of integers or a
            boolean mask of length number of pixels
        color: a matplotlib conform color
            the color for the pixel highlighting
        linewidth: float
            linewidth of the highlighting in points
        alpha: 0 <= alpha <= 1
            The transparency
        """

        lw_array = np.zeros_like(self.image)
        lw_array[pixels] = linewidth
        pixel_highlighting = copy(self.pixels)
        pixel_highlighting.set_facecolor('none')
        pixel_highlighting.set_linewidth(lw_array)
        pixel_highlighting.set_alpha(alpha)
        pixel_highlighting.set_edgecolor(color)
        self.ax.add_collection(pixel_highlighting)
        return pixel_highlighting

    def annotate_tm_edge_label(self):
        """
        Annotate each of the TMs on the top and bottom of the camera
        """
        if self._mapping is not None:
            kw = dict(fontsize=6, color='black', ha='center')
            m = self._mapping
            pix_size = self._mapping.metadata['size']
            f_tm_top = lambda g: m.ix[m.ix[g.index]['row'].idxmax(), 'slot']
            f_tm_bottom = lambda g: m.ix[m.ix[g.index]['row'].idxmin(), 'slot']
            tm_top = np.unique(m.groupby('col').agg(f_tm_top)['slot'])
            tm_bottom = np.unique(m.groupby('col').agg(f_tm_bottom)['slot'])
            for tm in tm_top:
                df = m.loc[m['slot'] == tm]
                ypix = df['ypix'].max() + pix_size * 0.7
                xpix = df['xpix'].mean()
                tm_txt = "TM{:02d}".format(tm)
                self.ax.text(xpix, ypix, tm_txt, va='bottom', **kw)
            for tm in tm_bottom:
                df = m.loc[m['slot'] == tm]
                ypix = df['ypix'].min() - pix_size * 0.7
                xpix = df['xpix'].mean()
                tm_txt = "TM{:02d}".format(tm)
                self.ax.text(xpix, ypix, tm_txt, va='top', **kw)
        else:
            print("Cannot annotate, no mapping attached to class")

    def annotate_led_flasher(self):
        """
        Annotate each of the LED flashers in the four corners of the camera
        """
        if self._mapping is not None:
            pix_size = self._mapping.metadata['size']
            axl = self._mapping.metadata['fOTUpX_l']
            ayl = self._mapping.metadata['fOTUpY_l'] + 2 * pix_size

            dxl = [1, -1, 1, -1]
            dyl = [1, 1, -1, -1]
            for i, (dx, dy) in enumerate(zip(dxl, dyl)):
                x = axl * dx
                y = ayl * dy
                self.ax.add_patch(Circle((x, y), radius=0.01, color='red'))
                self.ax.text(x, y, f"{i}", fontsize=7, color='white',
                             ha='center', va='center')
        else:
            print("Cannot annotate, no mapping attached to class")

    @classmethod
    def from_mapping(cls, mapping, **kwargs):
        """
        Generate the class from a CHECLabPy mapping dataframe

        Parameters
        ----------
        mapping : `pandas.DataFrame`
            The mapping for the pixels stored in a pandas DataFrame. Can be
            obtained from either of these options:

            CHECLabPy.io.TIOReader.mapping
            CHECLabPy.io.ReaderR0.mapping
            CHECLabPy.io.ReaderR1.mapping
            CHECLabPy.io.DL1Reader.mapping
            CHECLabPy.utils.mapping.get_clp_mapping_from_tc_mapping
        kwargs
            Arguments passed to `CHECLabPy.plottong.setup.Plotter`

        Returns
        -------
        `CameraImage`

        """
        xpix = mapping['xpix'].values
        ypix = mapping['ypix'].values
        size = mapping.metadata['size']
        image = cls(xpix, ypix, size, **kwargs)
        image._mapping = mapping
        return image

    @classmethod
    def from_tc_mapping(cls, tc_mapping, **kwargs):
        """
        Generate the class using the TargetCalib Mapping Class
        Parameters
        ----------
        tc_mapping : `target_calib.Mapping`
        kwargs
            Arguments passed to `CHECLabPy.plottong.setup.Plotter`

        Returns
        -------
        `CameraImage`

        """
        mapping = get_clp_mapping_from_tc_mapping(tc_mapping)
        return cls.from_mapping(mapping, **kwargs)

    @classmethod
    def from_camera_version(cls, camera_version, single=False, **kwargs):
        """
        Generate the class using the camera version (required TargetCalib)

        Parameters
        ----------
        camera_version : str
            Version of the camera (e.g. "1.0.1" corresponds to CHEC-S)
        single : bool
            Designate if it is just a single module you wish to plot
        kwargs
            Arguments passed to `CHECLabPy.plottong.setup.Plotter`

        Returns
        -------
        `CameraImage`

        """
        from target_calib import CameraConfiguration
        config = CameraConfiguration(camera_version)
        tc_mapping = config.GetMapping(single)
        return cls.from_tc_mapping(tc_mapping, **kwargs)

class CameraImageImshow(Plotter):
    def __init__(self, row, col, n_rows, n_cols, **kwargs):
        """
        Create a camera-image plot using imshow (essentially a 2D histogram,
        therefore missing module gaps)

        Parameters
        ----------
        row : ndarray
            The row for each pixel/superpixel/TM
        col : ndarray
            The row for each pixel/superpixel/TM
        n_rows : int
            The number of rows of pixel/superpixel/TM on the camera
        n_cols : int
            The number of rows of pixel/superpixel/TM on the camera
        kwargs
            Arguments passed to `CHECLabPy.plottong.setup.Plotter`
        """
        super().__init__(**kwargs)

        self._image = None
        self._mapping = None
        self.colorbar = None
        self.autoscale = True

        self.row = row
        self.col = col
        self.n_rows = n_rows
        self.n_cols = n_cols
        self.n_pixels = self.row.size

        assert self.row.size == self.col.size

        data = np.ma.zeros((self.n_rows, self.n_cols))
        self.mask(data)
        self.camera = self.ax.imshow(data, origin='lower')
        self.ax.axis('off')

    @property
    def image(self):
        return self._image

    @image.setter
    def image(self, val):
        assert val.size == self.n_pixels

        self._image = val
        data = np.ma.zeros((self.n_rows, self.n_cols))
        self.mask(data)
        data[self.row, self.col] = val
        self.camera.set_data(data)
        if self.autoscale:
            self.camera.autoscale()

    def add_colorbar(self, label=''):
        self.colorbar = self.fig.colorbar(self.camera, label=label)

    def set_limits_minmax(self, zmin, zmax):
        """
        Set the color scale limits from min to max
        """
        self.camera.set_clim(zmin, zmax)
        self.autoscale = False

    def reset_limits(self):
        """
        Reset to auto color scale limits
        """
        self.autoscale = True
        self.camera.autoscale()

    def mask(self, data):
        if (self.n_rows == 48) & (self.n_cols == 48):
            data[0:8, 40:48] = np.ma.masked
            data[0:8, 0:8] = np.ma.masked
            data[40:48, 0:8] = np.ma.masked
            data[40:48, 40:48] = np.ma.masked
        elif (self.n_rows == 24) & (self.n_cols == 24):
            data[0:4, 20:24] = np.ma.masked
            data[0:4, 0:4] = np.ma.masked
            data[20:24, 0:4] = np.ma.masked
            data[20:24, 20:24] = np.ma.masked
        elif (self.n_rows == 6) & (self.n_cols == 6):
            data[0:1, 5:6] = np.ma.masked
            data[0:1, 0:1] = np.ma.masked
            data[5:6, 0:1] = np.ma.masked
            data[5:6, 5:6] = np.ma.masked

    def annotate_on_telescope_up(self):
        """
        Add an arrow indicating where "ON-Telescope-UP" is
        """
        if self._mapping is not None:
            axl = self._mapping.metadata['fOTUpCol_l']
            ayl = self._mapping.metadata['fOTUpRow_l']
            adx = self._mapping.metadata['fOTUpCol_u'] - axl
            ady = self._mapping.metadata['fOTUpRow_u'] - ayl
            text = "ON-Telescope UP"
            self.ax.arrow(axl, ayl, adx, ady, head_width=1, head_length=1,
                          fc='r', ec='r')
            self.ax.text(axl, ayl, text, fontsize=8, color='r',
                         ha='center', va='bottom')
        else:
            print("Cannot annotate, no mapping attached to class")

    def add_text_to_pixel(self, pixel, value, size=3, color='w', **kwargs):
        """
        Add a text label to a single pixel

        Parameters
        ----------
        pixel : int
        value : str or float
        size : int
            Font size
        color : str
            Color of the text
        kwargs
            Named arguments to pass to matplotlib.pyplot.text
        """
        pos_x = self.col[pixel]
        pos_y = self.row[pixel]
        self.ax.text(pos_x, pos_y, value,
                     fontsize=size, color=color, ha='center', **kwargs)

    def add_pixel_text(self, values, size=3, color='w', **kwargs):
        """
        Add a text label to each pixel

        Parameters
        ----------
        values : ndarray
        size : int
            Font size
        color : str
            Color of the text
        kwargs
            Named arguments to pass to matplotlib.pyplot.text
        """
        assert values.size == self.n_pixels
        for pixel in range(self.n_pixels):
            self.add_text_to_pixel(pixel, values[pixel], size, color, **kwargs)

    def annotate_tm_edge_label(self):
        """
        Annotate each of the TMs on the top and bottom of the camera
        """
        if self._mapping is not None:
            kw = dict(fontsize=6, color='black', ha='center')
            m = self._mapping
            n_rows = m.metadata['n_rows']
            f_tm_top = lambda g: m.ix[m.ix[g.index]['row'].idxmax(), 'slot']
            f_tm_bottom = lambda g: m.ix[m.ix[g.index]['row'].idxmin(), 'slot']
            tm_top = np.unique(m.groupby('col').agg(f_tm_top)['slot'])
            tm_bottom = np.unique(m.groupby('col').agg(f_tm_bottom)['slot'])
            for tm in tm_top:
                df = m.loc[m['slot'] == tm]
                row = df['row'].max() + n_rows/48 + 0.4
                col = df['col'].mean()
                tm_txt = "TM{:02d}".format(tm)
                self.ax.text(col, row, tm_txt, va='bottom', **kw)
            for tm in tm_bottom:
                df = m.loc[m['slot'] == tm]
                row = df['row'].min() - n_rows/48 - 0.4
                col = df['col'].mean()
                tm_txt = "TM{:02d}".format(tm)
                self.ax.text(col, row, tm_txt, va='top', **kw)
        else:
            print("Cannot annotate, no mapping attached to class")

    @classmethod
    def from_mapping(cls, mapping, **kwargs):
        """
        Generate the class using a CHECLabPy mapping dataframe

        Parameters
        ----------
        mapping : `pandas.DataFrame`
            The mapping for the pixels stored in a pandas DataFrame. Can be
            obtained from either of these options:

            CHECLabPy.io.TIOReader.mapping
            CHECLabPy.io.ReaderR0.mapping
            CHECLabPy.io.ReaderR1.mapping
            CHECLabPy.io.DL1Reader.mapping
            CHECLabPy.utils.mapping.get_clp_mapping_from_tc_mapping
        kwargs
            Arguments passed to `CHECLabPy.plottong.setup.Plotter`

        Returns
        -------
        `CameraImage`

        """
        row = mapping['row'].values
        col = mapping['col'].values
        n_rows = mapping.metadata['n_rows']
        n_cols = mapping.metadata['n_columns']
        image = cls(row, col, n_rows, n_cols, **kwargs)
        image._mapping = mapping
        return image

    @classmethod
    def from_tc_mapping(cls, tc_mapping, **kwargs):
        """
        Generate the class using the TargetCalib Mapping Class
        Parameters
        ----------
        tc_mapping : `target_calib.Mapping`
        kwargs
            Arguments passed to `CHECLabPy.plottong.setup.Plotter`

        Returns
        -------
        `CameraImage`

        """
        mapping = get_clp_mapping_from_tc_mapping(tc_mapping)
        return cls.from_mapping(mapping, **kwargs)

    @classmethod
    def from_camera_version(cls, camera_version, single=False, **kwargs):
        """
        Generate the class using the camera version (required TargetCalib)

        Parameters
        ----------
        camera_version : str
            Version of the camera (e.g. "1.0.1" corresponds to CHEC-S)
        single : bool
            Designate if it is just a single module you wish to plot
        kwargs
            Arguments passed to `CHECLabPy.plottong.setup.Plotter`

        Returns
        -------
        `CameraImage`

        """
        from target_calib import CameraConfiguration
        config = CameraConfiguration(camera_version)
        tc_mapping = config.GetMapping(single)
        return cls.from_tc_mapping(tc_mapping, **kwargs)

max_stars_repo_name: ConteFrancesco/CHECLabPy

max_issues_repo_name: watsonjj/CHECLabPy

max_forks_repo_name: watsonjj/CHECLabPy


text:
import utils
from sagan_models import Generator, Discriminator
import paddle
import numpy as np
import random
from paddle.vision.transforms import functional as F
import time
import datetime
from PIL import Image
import os

class Tester(object):

    def __init__(self, config):

        # Config
        self.config = config

        # Create directories if not exist
        utils.make_folder(self.config.test_path)

        # Make dataloader
        self.dataloader, self.num_of_classes = utils.make_test_dataloader(batch_size=self.config.batch_size_test,
                                                                          dataset_type=self.config.dataset,
                                                                          data_path=self.config.test_data_path,
                                                                          shuffle=True,
                                                                          drop_last=self.config.drop_last,
                                                                          dataloader_args=self.config.dataloader_args,
                                                                          resize=self.config.resize,
                                                                          imsize=self.config.imsize,
                                                                          centercrop=self.config.centercrop,
                                                                          centercrop_size=self.config.centercrop_size,
                                                                          normalize=self.config.normalize,
                                                                          num_workers=self.config.num_workers,
                                                                          )

        # Data iterator
        self.data_iter = iter(self.dataloader)
        # Build G and D
        self.build_models()

    def build_models(self):
        self.G = Generator(self.config.z_dim, self.config.g_conv_dim, self.num_of_classes)
        self.D = Discriminator(self.config.d_conv_dim, self.num_of_classes)

        # Start with pretrained model (if it exists)
        if self.config.pretrained_model != '':
            utils.load_test_pretrained_model(self)

    def test(self):
        # For BatchNorm
        self.G.eval()

        # Init
        start_time = time.time()

        # Start testing
        for i in range(50000):
            print("id =", i)

            # Get real images & real labels (only need real labels)
            real_images, real_labels = self.get_real_samples()
            print("real_labels:", real_labels)

            # Create random noise
            z = paddle.randn([self.config.batch_size_test, self.config.z_dim])

            # Generate fake images for same real labels
            fake_images = self.G(z, real_labels)

            # Print out log info
            curr_time = time.time()
            curr_time_str = datetime.datetime.fromtimestamp(curr_time).strftime('%Y-%m-%d %H:%M:%S')
            elapsed = str(datetime.timedelta(seconds=(curr_time - start_time)))
            log = ("[{}] : Elapsed [{}]\n".
                   format(curr_time_str, elapsed))
            print('\n' + log)

            # Sample images
            print("Saving image samples..")
            self.G.eval()
            self.save_sample(fake_images, os.path.join(self.config.test_path + '/' + 'fake_{:05d}.png'.format(i)))

    def get_real_samples(self):
        try:
            real_images, real_labels = next(self.data_iter)
        except:
            self.data_iter = iter(self.dataloader)
            real_images, real_labels = next(self.data_iter)

        return real_images, real_labels

    def save_sample(self, images, path):
        # images = paddle.nn.functional.pad(images, pad = [2, 2, 2, 2])

        b, c, h, w = images.shape
        results = np.zeros((1, 3, 1 * h, 1 * w))
        count = 0
        for i in range(1):
            for j in range(1):
                results[:, :, i * h:(i + 1) * h, j * w:(j + 1) * w] = images[count].unsqueeze(0)
                count += 1
        results = 255 * (results + 1) / 2
        result = np.array(results[0].transpose(1, 2, 0), dtype=np.uint8)

        save_result = Image.fromarray(result)
        save_result.save(path)

max_stars_repo_name: Atmosphere-art/Self-Attention-GAN

max_issues_repo_name: Atmosphere-art/Self-Attention-GAN

max_forks_repo_name: Atmosphere-art/Self-Attention-GAN


text:
"""Generalized Gauss-Laguerre quadrature rule."""
import numpy
from scipy.special import gamma
import chaospy

from .hypercube import hypercube_quadrature

def laguerre(order, alpha=0., physicist=False):
    r"""
    Generalized Gauss-Laguerre quadrature rule.

    Compute the sample points and weights for Gauss-Laguerre quadrature. The
    sample points are the roots of the nth degree Laguerre polynomial. These
    sample points and weights correctly integrate polynomials of degree
    :math:`2N-1` or less.

    Gaussian quadrature come in two variants: physicist and probabilist. For
    Gauss-Laguerre physicist means a weight function :math:`x^\alpha e^{-x}`
    and weights that sum to :math`\Gamma(\alpha+1)`, and probabilist means a
    weight function is :math:`x^\alpha e^{-x}` and sum to 1.

    Args:
        order (int):
            The quadrature order.
        alpha (float):
            Shape parameter. Defaults to non-generalized Laguerre if 0.
        physicist (bool):
            Use physicist weights instead of probabilist.

    Returns:
        abscissas (numpy.ndarray):
            The ``order+1`` quadrature points for where to evaluate the model
            function with.
        weights (numpy.ndarray):
            The quadrature weights associated with each abscissas.

    Examples:
        >>> abscissas, weights = chaospy.quadrature.laguerre(2)
        >>> abscissas
        array([[0.41577456, 2.29428036, 6.28994508]])
        >>> weights
        array([0.71109301, 0.27851773, 0.01038926])

    See also:
        :func:`chaospy.quadrature.gaussian`

    """
    order = int(order)
    coefficients = chaospy.construct_recurrence_coefficients(
        order=order, dist=chaospy.Gamma(alpha+1))
    [abscissas], [weights] = chaospy.coefficients_to_quadrature(coefficients)
    weights *= gamma(alpha+1) if physicist else 1
    return abscissas[numpy.newaxis], weights

reason: import numpy,from scipy


text:
import numpy as np

color2num = dict(
    gray=30,
    red=31,
    green=32,
    yellow=33,
    blue=34,
    magenta=35,
    cyan=36,
    white=37,
    crimson=38
)

def colorize(string, color, bold=False, highlight=False):
    """
    Colorize a string.

    This function was originally written by John Schulman.
    """
    attr = []
    num = color2num[color]
    if highlight: num += 10
    attr.append(str(num))
    if bold: attr.append('1')
    return f'\x1b[{";".join(attr)}m{string}\x1b[0m'

def pwc(*args, color='red', bold=False, highlight=False):
    """
    Print with color
    """
    if isinstance(args, (tuple, list)):
        for s in args:
            print(colorize(s, color, bold, highlight))
    else:
        print(colorize(args, color, bold, highlight))

def assert_colorize(cond, err_msg=''):
    assert cond, colorize(err_msg, 'red')

def display_var_info(tf_vars, name='trainable', prefix=''):
    pwc(f'{prefix}{name.title()} variables', color='yellow')
    nparams = 0
    for v in tf_vars:
        name = v.name
        if '/Adam' in name or 'beta1_power' in name or 'beta2_power' in name: continue
        v_params = int(np.prod(v.shape.as_list()))
        nparams += v_params
        if '/b:' in name or '/biases' in name: continue    # Wx+b, bias is not interesting to look at => count params, but not print
        pwc(f'{prefix}   {name}{" "*(100-len(name))} {v_params:d} params {v.shape}', color='yellow')

    pwc(f'{prefix}Total model parameters: {nparams*1e-6:0.4g} million', color='yellow')

    return nparams

def display_model_var_info(models):
    learnable_models = {}
    opts = {}
    nparams = 0
    for name, model in models.items():
        if 'target' in name or name in learnable_models or name in opts:
            pass # ignore variables in the target networks
        elif 'opt' in name:
            opts[name] = model
        else:
            learnable_models[name] = model

    pwc(f'Learnable models:', color='yellow')
    for name, model in learnable_models.items():
        nparams += display_var_info(
            model.trainable_variables, name=name, prefix='   ')
    pwc(f'Total learnable model parameters: {nparams*1e-6:0.4g} million', color='yellow')

max_stars_repo_name: xlnwel/grl


text:
"""
This module provides functions for constructing Theano and Opendeep variables.
"""
# standard imports
import logging
# third-party libraries
import numpy
import theano
import theano.tensor as T
from six import integer_types
# internal imports
from opendeep.utils.config import create_dictionary_like

log = logging.getLogger(__name__)

def function(*args, **kwargs):
    """
    A wrapper around theano.function that disables the on_unused_input error (sets `on_unused_input` = "warn").
    Almost no part of OpenDeep can assume that an unused input is an error, so
    the default from Theano is inappropriate for this project.

    See: http://deeplearning.net/software/theano/library/compile/function.html

    Parameters
    ----------
    *args
        Variable length argument list to theano.function.
    **kwargs
        Arbitrary keyword arguments to theano.function.

    Returns
    -------
    theano.function
        Compiled Theano function.
    """
    return theano.function(*args, on_unused_input='warn', **kwargs)

def grad(*args, **kwargs):
    """
    A wrapper around theano.gradient.grad that disable the disconnected_inputs
    error (sets `disconnected_inputs` = "warn"). Almost no part of OpenDeep can assume that a disconnected input
    is an error.

    See: http://deeplearning.net/software/theano/library/tensor/basic.html#theano.gradient.grad

    Parameters
    ----------
    *args
        Variable length argument list to theano.gradient.grad.
    **kwargs
        Arbitrary keyword arguments to theano.gradient.grad.

    Returns
    -------
    variable or list/tuple of Variables (matching `wrt`)
        Symbolic expression of gradient of `cost` with respect to each of the `wrt` terms.
        If an element of `wrt` is not differentiable with respect to the output, then a zero variable is returned.
        It returns an object of same type as `wrt`: a list/tuple or Variable in all cases.
    """
    return theano.gradient.grad(*args, disconnected_inputs='warn', **kwargs)

def sharedX(value, name=None, borrow=False, dtype=theano.config.floatX):
    """
    Transform value into a theano shared variable of type `dtype` or theano.config.floatX

    For Theano shared variables, see: http://deeplearning.net/software/theano/library/compile/shared.html

    Parameters
    ----------
    value : number, array, vector, matrix, or tensor
        The input value to create into a Theano shared variable.
    name : string, optional
        The name for this shared variable.
    borrow : bool, optional
        The boolean `borrow` value to use in the Theano `shared` function.
    dtype : string, optional
        The `dtype` to use during conversion. Defaults to theano.config.floatX.

    Returns
    -------
    SharedVariable
        The Theano shared variable of the input `value`.
    """
    return theano.shared(theano._asarray(value, dtype=dtype),
                         name=name,
                         borrow=borrow)

def dataset_shared(dataset, name=None, borrow=False, dtype=theano.config.floatX):
    """
    Transform input `dataset` into a Theano shared variable of type `dtype`.

    .. todo:: Currently acts as a wrapper for `sharedX`. This is used for datasets,
        so we might want to use theano.tensor._shared instead of theano.shared for GPU optimizations.

    Parameters
    ----------
    dataset : number, array, vector, matrix, or tensor
        The input dataset to create into a Theano shared variable.
    name : string, optional
        The name for this shared variable.
    borrow : bool, optional
        The boolean `borrow` value to use in the Theano `shared` function.
    dtype : string, optional
        The `dtype` to use during conversion. Defaults to theano.config.floatX.

    Returns
    -------
    SharedVariable
        The Theano shared variable of the input `dataset`.
    """
    # try:
    #     return sharedX(value=dataset, name=name, borrow=borrow, dtype=dtype)
    # except MemoryError:
    #     warnings.warn("Dataset was too big to fit in single shared variable, returning a tensor._shared instead...")
    #     return theano.tensor._shared(value=dataset, name=name, borrow=borrow)
    return theano.tensor._shared(value=dataset, name=name, borrow=borrow)

def as_floatX(variable):
    """
    Casts a given variable into dtype `theano.config.floatX`. Numpy ndarrays will
    remain numpy ndarrays, python floats will become 0-D ndarrays and
    all other types will be treated as theano tensors.

    Parameters
    ----------
    variable: int, float, numpy array, or tensor
        The input to convert to type `theano.config.floatX`.

    Returns
    -------
    numpy array or tensor
        The input `variable` casted as type `theano.config.floatX`.
    """
    if isinstance(variable, (integer_types, float, numpy.number, numpy.ndarray)):
        return numpy.cast[theano.config.floatX](variable)

    return theano.tensor.cast(variable, theano.config.floatX)

def constantX(value):
    """
    Returns a constant of value `value` with type theano.config.floatX.

    Parameters
    ----------
    value : array_like
        The input value to create a Theano constant.

    Returns
    -------
    Constant
        A theano.tensor.constant of the input `value`.
    """
    return theano.tensor.constant(
        numpy.asarray(value, dtype=theano.config.floatX)
    )

def init_from_config(class_type, config):
    """
    Takes a class type and a configuration, and instantiates the class with the parameters in config.

    Parameters
    ----------
    class_type : Class
        The class type to create a new instance of.
    config : a json, yaml, or dictionary-like object
        The configuration to instantiate a new class of `class_type` with.

    Returns
    -------
    Class
        A newly instantiated class with parameters from `config`.
    """
    config_dict = create_dictionary_like(config)
    if config_dict is None:
        config_dict = {}
    return class_type(**config_dict)

max_stars_repo_name: vitruvianscience/OpenDeep

max_issues_repo_name: afcarl/OpenDeep

max_forks_repo_name: afcarl/OpenDeep

reason: import numpy,import theano


text:
import argparse
import random
from onto import Ontology
import numpy as np
import os
import json
from scipy import stats

def normalize(ont, hpid_filename, column=0):
    concepts = [c.strip().split()[column].replace("_",":") for c in open(hpid_filename).readlines() if c.strip()!=""]
    filtered = [ont.real_id[c] for c in concepts if c in ont.real_id]
    # and x.replace("_",":").strip()!="HP:0003220" and x.replace("_",":").strip()!="HP:0001263"and x.replace("_",":").strip()!="HP:0001999"]
    #raw = [ont.real_id[x.replace("_",":").strip()] for x in open(hpid_filename).readlines()]
    return set([c for c in filtered if c in ont.concepts])

def get_all_ancestors(ont, hit_list):
    return set([ont.concepts[x] for hit in hit_list for x in ont.ancestor_weight[ont.concept2id[hit]]])

def get_tp_fp(positives, real_positives):
    tp = len(positives & real_positives)
    fp = len(positives) - tp
    return tp, fp

def get_fmeasure(precision, recall):
    return 2.0*precision*recall/(precision+recall) if (precision+recall)!=0 else 0.0

def get_micro_stats(matrix):
    tp = matrix['tp']
    fp = matrix['fp']
    rp = matrix['rp']
    if np.sum(tp)+np.sum(fp) == 0:
        precision = 1.0
    else:
        precision = np.sum(tp)/(np.sum(tp)+np.sum(fp))
    if np.sum(rp) == 0:
        recall = 1.0
    else:
        recall = np.sum(tp)/np.sum(rp)
    return {"precision":precision, "recall":recall,
            "fmeasure":get_fmeasure(precision, recall)}

def get_macro_stats(matrix):
    tp = matrix['tp']
    fp = matrix['fp']
    rp = matrix['rp']
    precision = np.mean(np.where(tp+fp>0, tp/(tp+fp), 1.0))
    #precision = np.mean(np.where(tp+fp>0, tp/(tp+fp), 0.0))
    recall = np.mean(np.where(rp>0, tp/rp, 1.0))
    return {"precision":precision, "recall":recall,
            "fmeasure":get_fmeasure(precision, recall)}

def get_extended_stats(matrix):
    tp = matrix['tp']
    fp = matrix['fp']
    rp = matrix['rp']
    tp_precision = matrix['tp_ont_precision']
    tp_recall = matrix['tp_ont_recall']

    precision = np.mean(np.where(tp+fp>0, tp_precision/(tp+fp), 1.0))
    recall = np.mean(np.where(rp>0, tp_recall/rp, 1.0))
    return {"precision":precision, "recall":recall,
            "fmeasure":get_fmeasure(precision, recall)}

def print_results(results, is_mimic=False):
    res_print = []
    styles = ["micro", "macro"]
    if not is_mimic:
        styles.append("ont")
    for style in styles:
        for acc_type in ["precision", "recall", "fmeasure"]:
            res_print.append(results[style][acc_type])
    if not is_mimic:
        res_print.append(results['jaccard'])

    res_print = [x*100 for x in res_print]
    if is_mimic:
        print("%.1f\t%.1f\t%.1f\t%.1f\t%.1f\t%.1f" % tuple(res_print))
    else:
        print("%.1f\t%.1f\t%.1f\t%.1f\t%.1f\t%.1f\t%.1f\t%.1f\t%.1f\t%.1f" % tuple(res_print))

def get_confusion_matrix_mimic(label_dir, output_dir, file_list, ont, snomed2icd, column=0):
    true_positives = []
    false_positives = []
    real_positives = []
    total_calls = 0

    for filename in file_list:
        filename = filename.strip()

        file_real_positives = set([x.strip() for x in open(label_dir+"/"+filename).readlines() if x.strip() in snomed2icd.values()])

        file_positives = normalize(ont, output_dir+"/"+filename, column)
        file_positives = set([snomed2icd[x] for x in file_positives if x in snomed2icd])

        total_calls += len(file_positives)

        tp, fp = get_tp_fp(file_positives, file_real_positives)
        true_positives.append(tp)
        false_positives.append(fp)
        real_positives.append(len(file_real_positives))

    tp = np.array(true_positives)
    fp = np.array(false_positives)
    rp = np.array(real_positives)

    matrix = {
            'tp':np.array(true_positives),
            'fp':np.array(false_positives),
            'rp':np.array(real_positives),
            'total_calls': total_calls
            }
    return matrix

def get_confusion_matrix(label_dir, output_dir, file_list, ont, column=0):
    true_positives = []
    false_positives = []
    real_positives = []

    tp_ont_recall_list = []
    tp_ont_precision_list = []

    jaccard = []

    total_calls = 0

    for filename in file_list:
        filename = filename.strip()

        file_real_positives = normalize(ont, label_dir+"/"+filename)
        file_real_positives_ont = get_all_ancestors(ont, file_real_positives)

        file_positives = normalize(ont, output_dir+"/"+filename, column)
        total_calls += len(file_positives)
        file_positives_ont = get_all_ancestors(ont, file_positives)

        tp, fp = get_tp_fp(file_positives, file_real_positives)
        true_positives.append(tp)
        false_positives.append(fp)
        real_positives.append(len(file_real_positives))

        tp_ont_recall, _ = get_tp_fp(file_positives_ont, file_real_positives)
        tp_ont_precision, _ = get_tp_fp(file_positives, file_real_positives_ont)
        tp_ont_recall_list.append(tp_ont_recall)
        tp_ont_precision_list.append(tp_ont_precision)

        if len(file_real_positives)==0:
            jaccard.append(1.0)
        else:
            jaccard.append(
                    len(file_real_positives_ont & file_positives_ont)/
                    len(file_real_positives_ont | file_positives_ont))

    tp = np.array(true_positives)
    fp = np.array(false_positives)
    rp = np.array(real_positives)

    matrix = {
            'tp':np.array(true_positives),
            'fp':np.array(false_positives),
            'rp':np.array(real_positives),
            'tp_ont_precision':np.array(tp_ont_precision_list),
            'tp_ont_recall':np.array(tp_ont_recall_list),
            'jaccard': np.mean(jaccard),
            'total_calls': total_calls
            }
    return matrix

def eval(label_dir, output_dir, file_list, ont, column=0):
    matrix = get_confusion_matrix(label_dir, output_dir, file_list, ont, column)
    ret = { "ont": get_extended_stats(matrix),
            "macro": get_macro_stats(matrix),
            "micro": get_micro_stats(matrix),
            "jaccard":matrix['jaccard']}
    return ret

def eval_mimic(label_dir, output_dir, file_list, ont, snomed2icd, column=0):
    matrix = get_confusion_matrix_mimic(label_dir, output_dir, file_list, ont, snomed2icd, column)
    ret = { "macro": get_macro_stats(matrix),
            "micro": get_micro_stats(matrix)}
    return ret

def main():
    parser = argparse.ArgumentParser(description='Hello!')
    parser.add_argument('label_dir', help="Path to the directory where the input text files are located")
    parser.add_argument('output_dir', help="Path to the directory where the output files will be stored")
    parser.add_argument('--obofile', help="address to the ontology .obo file")
    parser.add_argument('--snomed2icd', help="address to the ontology .obo file")
    parser.add_argument('--oboroot', help="the concept in the ontology to be used as root (only this concept and its descendants will be used)")
    parser.add_argument('--file_list', help="Path to the directory where the output files will be stored")
    parser.add_argument('--comp_dir', help="Path to the directory where the output files will be stored")
    parser.add_argument('--no_error', action="store_true")
    parser.add_argument('--eval_mimic', action="store_true")
    parser.add_argument('--output_column', type=int, help="", default=0)
    args = parser.parse_args()
    if args.no_error:
        np.seterr(divide='ignore', invalid='ignore')

    if args.snomed2icd != None:
        with open(args.snomed2icd, 'r') as fp:
            snomed2icd = json.load(fp)

    file_list = os.listdir(args.label_dir)
    if args.file_list != None:
        file_list = [x.strip() for x in open(args.file_list).readlines()]

    ont = Ontology(args.obofile, args.oboroot)

    if args.eval_mimic:
        results = eval_mimic(args.label_dir, args.output_dir, file_list, ont, snomed2icd, column=args.output_column)
    else:
        results = eval(args.label_dir, args.output_dir, file_list, ont, column=args.output_column)
    print_results(results, args.eval_mimic)

if __name__ == "__main__":
	main()


text:
# coding: utf-8

# In[1]:

import sys
from utils import *
import numpy as np
import argparse
import os
import pickle
import cv2
from tqdm import tqdm
from loader import TestLoader
sys.path.append('../')
from train.model import P_Net,R_Net,O_Net
import train.config as config
from detection.detector import Detector
from detection.fcn_detector import FcnDetector
from detection.MtcnnDetector import MtcnnDetector

# In[3]:

def main(args):
    '''通过PNet或RNet生成下一个网络的输入'''
    size=args.input_size
    batch_size=config.batches
    min_face_size=config.min_face
    stride=config.stride
    thresh=config.thresh
    #模型地址
    model_path=['../model/PNet/','../model/RNet/','../model/ONet']
    if size==12:
        net='PNet'
        save_size=24
    elif size==24:
        net='RNet'
        save_size=48
    #图片数据地址
    base_dir='../data/WIDER_train/'
    #处理后的图片存放地址
    data_dir='../data/%d'%(save_size)
    neg_dir=os.path.join(data_dir,'negative')
    pos_dir=os.path.join(data_dir,'positive')
    part_dir=os.path.join(data_dir,'part')
    for dir_path in [neg_dir,pos_dir,part_dir]:
        if not os.path.exists(dir_path):
            os.makedirs(dir_path)
    detectors=[None,None,None]
    PNet=FcnDetector(P_Net,model_path[0])
    detectors[0]=PNet
    if net=='RNet':
        RNet=Detector(R_Net,24,batch_size[1],model_path[1])
        detectors[1]=RNet
    basedir='../data/'
    filename='../data/wider_face_train_bbx_gt.txt'
    #读取文件的image和box对应函数在utils中
    data=read_annotation(base_dir,filename)
    mtcnn_detector=MtcnnDetector(detectors,min_face_size=min_face_size,
                                stride=stride,threshold=thresh)
    save_path=data_dir
    save_file=os.path.join(save_path,'detections.pkl')
    if not os.path.exists(save_file):
        #将data制作成迭代器
        print('载入数据')
        test_data=TestLoader(data['images'])
        detectors,_=mtcnn_detector.detect_face(test_data)
        print('完成识别')

        with open(save_file,'wb') as f:
            pickle.dump(detectors,f,1)
    print('开始生成图像')
    save_hard_example(save_size,data,neg_dir,pos_dir,part_dir,save_path)

# In[2]:

def save_hard_example(save_size, data,neg_dir,pos_dir,part_dir,save_path):
    '''将网络识别的box用来裁剪原图像作为下一个网络的输入'''

    im_idx_list = data['images']

    gt_boxes_list = data['bboxes']
    num_of_images = len(im_idx_list)

    # save files
    neg_label_file = "../data/%d/neg_%d.txt" % (save_size, save_size)
    neg_file = open(neg_label_file, 'w')

    pos_label_file = "../data/%d/pos_%d.txt" % (save_size, save_size)
    pos_file = open(pos_label_file, 'w')

    part_label_file = "../data/%d/part_%d.txt" % (save_size, save_size)
    part_file = open(part_label_file, 'w')
    #read detect result
    det_boxes = pickle.load(open(os.path.join(save_path, 'detections.pkl'), 'rb'))
    # print(len(det_boxes), num_of_images)

    assert len(det_boxes) == num_of_images, "弄错了"

    n_idx = 0
    p_idx = 0
    d_idx = 0
    image_done = 0

    for im_idx, dets, gts in tqdm(zip(im_idx_list, det_boxes, gt_boxes_list)):
        gts = np.array(gts, dtype=np.float32).reshape(-1, 4)
        image_done += 1

        if dets.shape[0] == 0:
            continue
        img = cv2.imread(im_idx)
        #转换成正方形
        dets = convert_to_square(dets)
        dets[:, 0:4] = np.round(dets[:, 0:4])
        neg_num = 0
        for box in dets:
            x_left, y_top, x_right, y_bottom, _ = box.astype(int)
            width = x_right - x_left + 1
            height = y_bottom - y_top + 1

            # 除去过小的
            if width < 20 or x_left < 0 or y_top < 0 or x_right > img.shape[1] - 1 or y_bottom > img.shape[0] - 1:
                continue

            Iou = IOU(box, gts)
            cropped_im = img[y_top:y_bottom + 1, x_left:x_right + 1, :]
            resized_im = cv2.resize(cropped_im, (save_size, save_size),
                                    interpolation=cv2.INTER_LINEAR)

            #划分种类
            if np.max(Iou) < 0.3 and neg_num < 60:

                save_file = os.path.join(neg_dir, "%s.jpg" % n_idx)

                neg_file.write(save_file + ' 0\n')
                cv2.imwrite(save_file, resized_im)
                n_idx += 1
                neg_num += 1
            else:

                idx = np.argmax(Iou)
                assigned_gt = gts[idx]
                x1, y1, x2, y2 = assigned_gt

                #偏移量
                offset_x1 = (x1 - x_left) / float(width)
                offset_y1 = (y1 - y_top) / float(height)
                offset_x2 = (x2 - x_right) / float(width)
                offset_y2 = (y2 - y_bottom) / float(height)

                # pos和part
                if np.max(Iou) >= 0.65:
                    save_file = os.path.join(pos_dir, "%s.jpg" % p_idx)
                    pos_file.write(save_file + ' 1 %.2f %.2f %.2f %.2f\n' % (
                        offset_x1, offset_y1, offset_x2, offset_y2))
                    cv2.imwrite(save_file, resized_im)
                    p_idx += 1

                elif np.max(Iou) >= 0.4:
                    save_file = os.path.join(part_dir, "%s.jpg" % d_idx)
                    part_file.write(save_file + ' -1 %.2f %.2f %.2f %.2f\n' % (
                        offset_x1, offset_y1, offset_x2, offset_y2))
                    cv2.imwrite(save_file, resized_im)
                    d_idx += 1
    neg_file.close()
    part_file.close()
    pos_file.close()

# In[4]:

def parse_arguments(argv):

    parser = argparse.ArgumentParser()

    parser.add_argument('input_size', type=int,
                        help='The input size for specific net')

    return parser.parse_args(argv)

if __name__ == '__main__':
    main(parse_arguments(sys.argv[1:]))