text:
import collections
import numpy as np

from menpo.shape import TriMesh
from menpo.image import Image

def check_param(n, types, param_name, param):
    r"""
    Utility method for validating a given parameter. Since parameters may be
    a single option applied to all items in the list, or a list of different
    options per item, we validate this. If a single argument is passed, it
    is extended into a list.

    Parameters
    ----------
    n : `int`
        Number of items expected.
    types : `type`
        Python `type` or list of `types` that the parameters are expected to be.
    param_name : `str`
        The name of parameter for printing error messages.
    param : `object`
        An object or list of objects that represent a parameter for a
        function.

    Returns
    -------
    param_list : `list` of type `object`
        This is a validated list of length ``n``, one parameter for each
        element in the implicit list of other parameters.
    """
    error_msg = "{0} must be a {1} or a list/tuple of " "{1} with length {2}".format(
        param_name, types, n
    )

    # Could be a single value - or we have an error
    if isinstance(param, types):
        return [param] * n
    elif not isinstance(param, collections.Iterable):
        raise ValueError(error_msg)

    # Must be an iterable object
    len_param = len(param)
    isinstance_all_in_param = all(isinstance(p, types) for p in param)

    if len_param == 1 and isinstance_all_in_param:
        return list(param) * n
    elif len_param == n and isinstance_all_in_param:
        return list(param)
    else:
        raise ValueError(error_msg)

def _rasterize_matplotlib(
    image,
    pclouds,
    render_lines=True,
    line_style="-",
    line_colour="b",
    line_width=1,
    render_markers=True,
    marker_style="o",
    marker_size=1,
    marker_face_colour="b",
    marker_edge_colour="b",
    marker_edge_width=1,
):
    import matplotlib.pyplot as plt

    # TODO: Since upgrading to Matplotlib 2.0 it seems that the rendering acts
    #       strangely with respect to pixel accurate rendering. The lines
    #       appear to be rendering too short - but seem correct when we reset
    #       to the classic style - so obviously there is some matplotlib setting
    #       we could just manually set to get the correct behaviour but it is
    #       a chore to find so I'm pragmatically maintaining the old behaviour
    #       here with this contextmanager.
    with plt.style.context("classic"):
        # Convert image shape into 100 DPI inches
        # This makes sure we maintain the original image size
        image_shape = np.array(image.shape)[::-1] / 100.0
        f = plt.figure(figsize=image_shape, frameon=False, dpi=100)

        image.view(figure_id=f.number, figure_size=image_shape)
        for k, p in enumerate(pclouds):
            p.view(
                figure_id=f.number,
                render_axes=False,
                figure_size=image_shape,
                render_lines=render_lines[k],
                line_style=line_style[k],
                line_colour=line_colour[k],
                line_width=line_width[k],
                render_markers=render_markers[k],
                marker_style=marker_style[k],
                marker_size=marker_size[k],
                marker_face_colour=marker_face_colour[k],
                marker_edge_colour=marker_edge_colour[k],
                marker_edge_width=marker_edge_width[k],
            )

        # Make sure the layout is tight so that the image is of the original size
        f.tight_layout(pad=0)
        # Get the pixels directly from the canvas buffer which is fast
        c_buffer, shape = f.canvas.print_to_buffer()
        # Turn buffer into numpy array and reshape to image
        shape = shape[::-1] + (-1,)
        pixels_buffer = np.frombuffer(c_buffer, dtype=np.uint8).reshape(shape)
        # Prevent matplotlib from rendering
        plt.close(f)
        # We have to apply the alpha channel to get the correct colour
        mask = pixels_buffer[..., 3:] * (1.0 / 255.0)
        pixels = (pixels_buffer[..., :3] * mask).astype(np.uint8)
        return Image.init_from_channels_at_back(pixels)

def _parse_colour(x):
    if isinstance(x, str):
        from matplotlib.colors import ColorConverter

        c = ColorConverter()
        x = c.to_rgb(x)

    if isinstance(x, (tuple, list)):
        # Assume we have a floating point rgb
        if any(a <= 1.0 for a in x):
            x = [int(a * 255) for a in x]
    return tuple(x)

def _rasterize_pillow(
    image,
    pclouds,
    render_lines=True,
    line_style="-",
    line_colour="b",
    line_width=1,
    render_markers=True,
    marker_style="o",
    marker_size=1,
    marker_face_colour="b",
    marker_edge_colour="b",
    marker_edge_width=1,
):
    from PIL import ImageDraw

    if image.n_channels not in {1, 3}:
        raise ValueError(
            "The Pillow rasterizer only supports grayscale or " "colour images"
        )
    if any(x != "-" for x in line_style):
        raise ValueError("The Pillow rasterizer only supports the '-' " "line style.")
    if any(x not in {"o", "s"} for x in marker_style):
        raise ValueError(
            "The Pillow rasterizer only supports the 'o' and 's' " "marker styles."
        )
    if any(x > 1 for x in marker_edge_width):
        raise ValueError(
            "The Pillow rasterizer only supports " "marker_edge_width of 1 or 0."
        )

    if image.n_channels == 1:
        # Make the image RGB
        image = image.extract_channels([0, 0, 0])

    pil_im = image.as_PILImage()
    draw = ImageDraw.Draw(pil_im)

    line_colour = [_parse_colour(x) for x in line_colour]
    marker_edge_colour = [_parse_colour(x) for x in marker_edge_colour]
    marker_face_colour = [_parse_colour(x) for x in marker_face_colour]

    for k in range(len(pclouds)):
        p = pclouds[k]
        if isinstance(p, TriMesh):
            pclouds[k] = p.as_pointgraph()

        points = p.points
        if (
            render_lines[k]
            and line_width[k] > 0
            and hasattr(p, "edges")
            and p.edges.size > 0
        ):
            edges = p.edges
            lines = zip(points[edges[:, 0], :], points[edges[:, 1], :])

            for l1, l2 in lines:
                draw.line(
                    [tuple(l1[::-1]), tuple(l2[::-1])],
                    fill=line_colour[k],
                    width=line_width[k],
                )

        if render_markers[k] and marker_size[k] > 0:
            draw_func = draw.ellipse if marker_style[k] == "o" else draw.rectangle
            outline = marker_edge_colour[k] if marker_edge_width[k] == 1 else None
            for p in points:
                y, x = p
                draw_func(
                    (
                        x - marker_size[k],
                        y - marker_size[k],
                        x + marker_size[k],
                        y + marker_size[k],
                    ),
                    fill=marker_face_colour[k],
                    outline=outline,
                )

    del draw

    pixels = np.asarray(pil_im)
    if image.n_channels == 3:
        return Image.init_from_channels_at_back(pixels)
    else:
        return Image(pixels)

_RASTERIZE_BACKENDS = {"matplotlib": _rasterize_matplotlib, "pillow": _rasterize_pillow}

def rasterize_landmarks_2d(
    image,
    group=None,
    render_lines=True,
    line_style="-",
    line_colour="b",
    line_width=1,
    render_markers=True,
    marker_style="o",
    marker_size=1,
    marker_face_colour="b",
    marker_edge_colour="b",
    marker_edge_width=1,
    backend="matplotlib",
):
    r"""
    This method provides the ability to rasterize 2D landmarks onto an image.
    The returned image has the specified landmark groups rasterized onto
    the image - which is useful for things like creating result examples
    or rendering videos with annotations.

    Since multiple landmark groups can be specified, all arguments can
    take lists of parameters that map to the provided groups list. Therefore,
    the parameters must be lists of the correct length or a single parameter
    to apply to every landmark group.

    Multiple backends are provided, all with different strengths. The 'pillow'
    backend is very fast, but not very flexible. The `matplotlib` backend
    should be feature compatible with other Menpo rendering methods, but
    is much slower due to the overhead of creating a figure to render
    into.

    Parameters
    ----------
    image : :map:`Image` or subclass
        The image to render onto.
    group : `str` or `list` of `str`, optional
        The landmark group key, or a list of keys.
    render_lines : `bool`, optional
        If ``True``, and the provided landmark group is a :map:`PointGraph`,
        the edges are rendered.
    line_style : `str`, optional
        The style of the edge line. Not all backends support this argument.
    line_colour : `str` or `tuple`, optional
        A Matplotlib style colour or a backend dependant colour.
    line_width : `int`, optional
        The width of the line to rasterize.
    render_markers : `bool`, optional
        If ``True``, render markers at the coordinates of each landmark.
    marker_style : `str`, optional
        A Matplotlib marker style. Not all backends support all marker styles.
    marker_size : `int`, optional
        The size of the marker - different backends use different scale
        spaces so consistent output may by difficult.
    marker_face_colour : `str`, optional
        A Matplotlib style colour or a backend dependant colour.
    marker_edge_colour : `str`, optional
        A Matplotlib style colour or a backend dependant colour.
    marker_edge_width : `int`, optional
        The width of the marker edge. Not all backends support this.
    backend : {'matplotlib', 'pillow'}, optional
        The backend to use.

    Returns
    -------
    rasterized_image : :map:`Image`
        The image with the landmarks rasterized directly into the pixels.
        The pixels of the image returned are of uint8 type.

    Raises
    ------
    ValueError
        Only 2D images are supported.
    ValueError
        Only RGB (3-channel) or Greyscale (1-channel) images are supported.
    """
    if image.n_channels != 1 and image.n_channels != 3:
        raise ValueError("Only RGB or Greyscale images can be rasterized")
    if image.n_dims != 2:
        raise ValueError("Only 2D images can be rasterized.")

    if backend in _RASTERIZE_BACKENDS:
        if isinstance(group, list):
            landmarks = [image.landmarks[g] for g in group]
        else:
            landmarks = [image.landmarks[group]]

        # Validate all the parameters for multiple landmark groups being
        # passed in
        n_pclouds = len(landmarks)
        render_lines = check_param(n_pclouds, bool, "render_lines", render_lines)
        line_style = check_param(n_pclouds, str, "line_style", line_style)
        line_colour = check_param(n_pclouds, (str, tuple), "line_colour", line_colour)
        line_width = check_param(n_pclouds, int, "line_width", line_width)
        render_markers = check_param(n_pclouds, bool, "render_markers", render_markers)
        marker_style = check_param(n_pclouds, str, "marker_style", marker_style)
        marker_size = check_param(n_pclouds, int, "marker_size", marker_size)
        marker_face_colour = check_param(
            n_pclouds, (str, tuple), "marker_face_colour", marker_face_colour
        )
        marker_edge_colour = check_param(
            n_pclouds, (str, tuple), "marker_edge_colour", marker_edge_colour
        )
        marker_edge_width = check_param(
            n_pclouds, int, "marker_edge_width", marker_edge_width
        )

        return _RASTERIZE_BACKENDS[backend](
            image,
            landmarks,
            render_lines=render_lines,
            line_style=line_style,
            line_colour=line_colour,
            line_width=line_width,
            render_markers=render_markers,
            marker_style=marker_style,
            marker_size=marker_size,
            marker_face_colour=marker_face_colour,
            marker_edge_colour=marker_edge_colour,
            marker_edge_width=marker_edge_width,
        )
    else:
        raise ValueError("Unsupported backend: {}".format(backend))

max_stars_repo_name: apapaion/menpo

max_issues_repo_name: apapaion/menpo

max_forks_repo_name: apapaion/menpo


text:
import numpy as np
import matplotlib.pyplot as plt
import cv2

def ImageInActualSize(I):
    dpi = plt.rcParams['figure.dpi']
    H,W = I.shape[:2]
    figSize = W/float(dpi) , H/float(dpi)
    fig = plt.figure(figsize = figSize)
    ax = fig.add_axes([0,0,1,1])
    ax.axis('off')
    ax.imshow(I,cmap='gray')
    plt.show

def resize(scale, I): # Using opencv function
    fx, fy = scale
    res = cv2.resize(I, fx=fx, fy=fy, dsize=None)
    return res

def bilinearInterpolate(r, c, I):
    lc = int(c)
    rc = lc + 1

    tr = int(r)
    br = tr + 1

    wl = rc - c
    wr = c - lc

    wt = br - r
    wb = r - tr

    if((tr>=0 and br<I.shape[0]) and (lc>=0 and rc<I.shape[1])):
        a = wl*I[tr, lc] + wr*I[tr, rc]
        b = wl*I[br, lc] + wr*I[br, rc]

        g = wt*a+wb*b

        return g
    else:
        return 0

def ScaleImageAlgo(scale, I):
    fx, fy = scale
    S = np.array([[fx, 0], [0, fy]])

    numRows = I.shape[0]
    numcols = I.shape[1]

    I2 = np.zeros((fx * numRows, fy * numcols), dtype="uint8")
    Tinv = np.linalg.inv(S)

    for i in range(I2.shape[0]):
        for j in range(I2.shape[1]):
            p_dash = np.array([i, j])
            p = Tinv.dot(p_dash)

            new_i, new_j = p
            if((new_i<0 or new_i>=numRows) or (new_j<0 or new_j>=numcols)):
                pass
            else:
                g = bilinearInterpolate(new_i, new_j, I)
                I2[i, j] = g
    return I2

# scale is a tuple consisting of width factor and height factor respectively.
# I is the numpy ndarray of image.

scale = (2, 2)
grayImage = r'../img/albert-einstein_gray.jpg'
colourImage = r'../img/tulips.jpg'
I_gray = cv2.imread(grayImage, cv2.IMREAD_GRAYSCALE)
I_BGR = cv2.imread(colourImage)

""" I_resized_gray = resize(scale, I_gray)
I_resized_BGR = resize(scale, I_BGR) """

I_resized_gray = ScaleImageAlgo(scale, I_gray)

I_col = I_BGR[:,:,::-1]
r = ScaleImageAlgo(scale, I_col[:,:,0])
g = ScaleImageAlgo(scale, I_col[:,:,1])
b = ScaleImageAlgo(scale, I_col[:,:,2])
C = np.zeros((r.shape[0], r.shape[1], 3), dtype="uint8")
C[:,:,0] = r
C[:,:,1] = g
C[:,:,2] = b

plt.figure(1)
plt.subplot(211)
plt.imshow(I_resized_gray, cmap = "gray")
plt.subplot(212)
plt.imshow(C)
plt.show()

ImageInActualSize(I_resized_gray)
plt.imsave(r'../img/resized_gray_using_algo.jpg', I_resized_gray, cmap = "gray")
ImageInActualSize(C)
plt.imsave(r'../img/resized_tulips_using_algo.jpg', C)

ext: py

reason: import numpy


text:
#!/usr/bin/env python3
# encoding: utf-8
# Code modified from https://github.com/Wangyixinxin/ACN
import glob
import os
import numpy as np
import nibabel as nib
import torch
from torch.utils.data import Dataset, DataLoader
import random

class Brats2018(Dataset):

    def __init__(self, patients_dir, crop_size, modes, train=True, normalization = True):
        self.patients_dir = patients_dir
        self.modes = modes
        self.train = train
        self.crop_size = crop_size
        self.normalization = normalization

    def __len__(self):
        return len(self.patients_dir)

    def __getitem__(self, index):
        patient_dir = self.patients_dir[index]
        volumes = []
        modes = list(self.modes) + ['seg']
        for mode in modes:
            patient_id = os.path.split(patient_dir)[-1]
            volume_path = os.path.join(patient_dir, patient_id + "_" + mode + '.nii')
            volume = nib.load(volume_path).get_data()
            if not mode == "seg" and self.normalization:
                volume = self.normlize(volume)  # [0, 1.0]
            volumes.append(volume)                  # [h, w, d]
        seg_volume = volumes[-1]
        volumes = volumes[:-1]
        volume, seg_volume = self.aug_sample(volumes, seg_volume)
        ed_volume = (seg_volume == 2) # peritumoral edema ED
        net_volume = (seg_volume == 1) # enhancing tumor core NET
        et_volume = (seg_volume == 4) # enhancing tumor ET
        bg_volume = (seg_volume == 0)

        seg_volume = [ed_volume, net_volume, et_volume, bg_volume]
        seg_volume = np.concatenate(seg_volume, axis=0).astype("float32")

        return (torch.tensor(volume.copy(), dtype=torch.float),
                torch.tensor(seg_volume.copy(), dtype=torch.float))

    def aug_sample(self, volumes, mask):
        """
            Args:
                volumes: list of array, [h, w, d]
                mask: array [h, w, d], segmentation volume
            Ret: x, y: [channel, h, w, d]

        """
        x = np.stack(volumes, axis=0)       # [N, H, W, D]
        y = np.expand_dims(mask, axis=0)    # [channel, h, w, d]

        if self.train:
            # crop volume
            x, y = self.random_crop(x, y)
            if random.random() < 0.5:
                x = np.flip(x, axis=1)
                y = np.flip(y, axis=1)
            if random.random() < 0.5:
                x = np.flip(x, axis=2)
                y = np.flip(y, axis=2)
            if random.random() < 0.5:
                x = np.flip(x, axis=3)
                y = np.flip(y, axis=3)
        else:
            x, y = self.center_crop(x, y)

        return x, y

    def random_crop(self, x, y):
        """
        Args:
            x: 4d array, [channel, h, w, d]
        """
        crop_size = self.crop_size
        height, width, depth = x.shape[-3:]
        sx = random.randint(0, height - crop_size[0] - 1)
        sy = random.randint(0, width - crop_size[1] - 1)
        sz = random.randint(0, depth - crop_size[2] - 1)
        crop_volume = x[:, sx:sx + crop_size[0], sy:sy + crop_size[1], sz:sz + crop_size[2]]
        crop_seg = y[:, sx:sx + crop_size[0], sy:sy + crop_size[1], sz:sz + crop_size[2]]

        return crop_volume, crop_seg

    def center_crop(self, x, y):
        crop_size = self.crop_size
        height, width, depth = x.shape[-3:]
        sx = (height - crop_size[0] - 1) // 2
        sy = (width - crop_size[1] - 1) // 2
        sz = (depth - crop_size[2] - 1) // 2
        crop_volume = x[:, sx:sx + crop_size[0], sy:sy + crop_size[1], sz:sz + crop_size[2]]
        crop_seg = y[:, sx:sx + crop_size[0], sy:sy + crop_size[1], sz:sz + crop_size[2]]

        return crop_volume, crop_seg

    def normlize(self, x):
        return (x - x.min()) / (x.max() - x.min())

    def normlize_brain(self, x, epsilon=1e-8):
        average        = x[np.nonzero(x)].mean()
        std            = x[np.nonzero(x)].std() + epsilon
        mask           = x>0
        sub_mean       = np.where(mask, x-average, x)
        x_normalized   = np.where(mask, sub_mean/std, x)
        return x_normalized

def split_dataset(data_root, test_p):
    patients_dir = glob.glob(os.path.join(data_root, "*GG", "Brats18*"))
    patients_dir.sort()
    N = int(len(patients_dir)*test_p)
    train_patients_list =  patients_dir[N:]
    val_patients_list   =  patients_dir[:N]

    return train_patients_list, val_patients_list

def make_data_loaders(config):
    train_list, val_list = split_dataset(config['path_to_data'], float(config['test_p']))
    crop_size = np.zeros((3))
    crop_size[0] = config['inputshape'][0]
    crop_size[1] = config['inputshape'][1]
    crop_size[2] = config['inputshape'][2]
    crop_size    = crop_size.astype(np.uint16)
    crop_size    = (160, 192, 128)
    train_ds = Brats2018(train_list, crop_size=crop_size, modes=config['modalities'], train=True)
    val_ds = Brats2018(val_list,     crop_size=crop_size, modes=config['modalities'], train=False)
    loaders = {}
    loaders['train'] = DataLoader(train_ds, batch_size=int(config['batch_size_tr']),
                                  num_workers=4,
                                  pin_memory=True,
                                  shuffle=True)
    loaders['eval'] = DataLoader(val_ds, batch_size=int(config['batch_size_va']),
                                  num_workers=4,
                                  pin_memory=True,
                                  shuffle=False)
    return loaders


text:
import numpy as np
import scipy.stats
from scipy.linalg import LinAlgError
from scipy.optimize import minimize
from multiprocessing import Pool
from functools import partial
import platform

from mogp_emulator.GaussianProcess import GaussianProcessBase, GaussianProcess
from mogp_emulator.GaussianProcessGPU import GaussianProcessGPU
from mogp_emulator.MultiOutputGP_GPU import MultiOutputGP_GPU
from mogp_emulator import LibGPGPU
from mogp_emulator.MultiOutputGP import MultiOutputGP

def fit_GP_MAP(*args, n_tries=15, theta0=None, method="L-BFGS-B",
               skip_failures=True, refit=False, **kwargs):
    """Fit one or more Gaussian Processes by attempting to minimize the
    negative log-posterior

    Fits the hyperparameters of one or more Gaussian Processes by
    attempting to minimize the negative log-posterior multiple times
    from a given starting location and using a particular minimization
    method. The best result found among all of the attempts is
    returned, unless all attempts to fit the parameters result in an
    error (see below).

    The arguments to the method can either be an existing
    ``GaussianProcess`` or ``MultiOutputGP`` instance, or a list of
    arguments to be passed to the ``__init__`` method of
    ``GaussianProcess`` or ``MultiOutputGP`` if more than one output
    is detected.  Keyword arguments for creating a new
    ``GaussianProcess`` or ``MultiOutputGP`` object can either be
    passed as part of the ``*args`` list or as keywords (if present in
    ``**kwargs``, they will be extracted and passed separately to the
    ``__init__`` method).

    If the method encounters an overflow (this can result because the
    parameter values stored are the logarithm of the actual
    hyperparameters to enforce positivity) or a linear algebra error
    (occurs when the covariance matrix cannot be inverted, even with
    additional noise added along the diagonal if adaptive noise was
    selected), the iteration is skipped. If all attempts to find
    optimal hyperparameters result in an error, then the emulator is
    skipped and the parameters are reset to ``None``. By default, a
    warning will be printed to the console if this occurs for
    ``MultiOutputGP`` fitting, while an error will be raise for
    ``GaussianProcess`` fitting. This behavior can be changed to raise
    an error for ``MultiOutputGP`` fitting by passing the kwarg
    ``skip_failures=False``. This default behavior is chosen because
    ``MultiOutputGP`` fitting is often done in a situation where human
    review of all fit emulators is not possible, so the fitting
    routine skips over failures and then flags those that failed for
    further review.

    For ``MultiOutputGP`` fitting, by default the routine assumes that
    only GPs that currently do not have hyperparameters fit need to be
    fit. This behavior is controlled by the ``refit`` keyword, which
    is ``False`` by default. To fit all emulators regardless of their
    current fitting status, pass ``refit=True``. The ``refit``
    argument has no effect on fitting of single ``GaussianProcess``
    objects -- standard ``GaussianProcess`` objects will be fit
    regardless of the current value of the hyperparameters.

    The ``theta0`` parameter is the point at which the first iteration
    will start. If more than one attempt is made, subsequent attempts
    will use random starting points. If you are fitting Multiple
    Outputs, then this argument can take any of the following forms:
    (1) None (random start points for all emulators, which are drawn
    from the prior distribution for each fit parameter), (2) a list of
    numpy arrays or ``NoneTypes`` with length ``n_emulators``, (3) a
    numpy array of shape ``(n_params,)`` or ``(n_emulators,
    n_params)`` which with either use the same start point for all
    emulators or the specified start point for all emulators. Note
    that if you us a numpy array, all emulators must have the same
    number of parameters, while using a list allows more flexibility.

    The user can specify the details of the minimization method, using
    any of the gradient-based optimizers available in
    ``scipy.optimize.minimize``. Any additional parameters beyond the
    method specification can be passed as keyword arguments.

    The function returns a fit ``GaussianProcess`` or
    ``MultiOutputGP`` instance, either the original one passed to the
    function, or the new one created from the included arguments.

    :param ``*args``: Either a single ``GaussianProcess`` or
                      ``MultiOutputGP`` instance, or arguments to be
                      passed to the ``__init__`` method when creating
                      a new ``GaussianProcess`` or ``MultiOutputGP``
                      instance.
    :param n_tries: Number of attempts to minimize the negative
                    log-posterior function.  Must be a positive
                    integer (optional, default is 15)
    :type n_tries: int
    :param theta0: Initial starting point for the first iteration. If
                   present, must be array-like with shape
                   ``(n_params,)`` based on the specific
                   ``GaussianProcess`` being fit. If a
                   ``MultiOutputGP`` is being fit it must be a list of
                   length ``n_emulators`` with each entry as either
                   ``None`` or a numpy array of shape ``(n_params,)``,
                   or a numpy array with shape ``(n_emulators,
                   n_params)`` (note that if the various emulators
                   have different numbers of parameters, the numpy
                   array option will not work).  If ``None`` is given,
                   then a random value is chosen. (Default is
                   ``None``)
    :type theta0: None or ndarray
    :param method: Minimization method to be used. Can be any
                   gradient-based optimization method available in
                   ``scipy.optimize.minimize``. (Default is
                   ``'L-BFGS-B'``)
    :type method: str
    :param skip_failures: Boolean controlling how to handle failures
                          in ``MultiOutputGP`` fitting. If set to
                          ``True``, emulator fits will fail silently
                          without raising an error and provide
                          information on the emulators that failed
                          and the end of fitting. If ``False``, any
                          failed fit will raise a ``RuntimeError``.
                          Has no effect on fitting a single
                          ``GaussianProcess``, which will always
                          raise an error. Optional, default is
                          ``True``.
    :type skip_failures: bool
    :param refit: Boolean indicating if previously fit emulators
                  for ``MultiOutputGP`` objects should be fit again.
                  Optional, default is ``False``. Has no effect on
                  ``GaussianProcess`` fitting, which will be fit
                  irrespective of the current hyperparameter values.
    :type refit: bool
    :param ``**kwargs``: Additional keyword arguments to be passed to
                         ``GaussianProcess.__init__``,
                         ``MultiOutputGP.__init__``, or the
                         minimization routine. Relevant parameters for
                         the GP classes are automatically split out
                         from those used in the minimization
                         function. See available parameters in the
                         corresponding functions for details.
    :returns: Fit GP or Multi-Output GP instance
    :rtype: GaussianProcess or MultiOutputGP or GaussianProcessGPU or MultiOutputGP_GPU

    """

    if len(args) == 1:
        gp = args[0]
        if isinstance(gp, MultiOutputGP):
            gp = _fit_MOGP_MAP(gp, n_tries, theta0, method, refit, **kwargs)
        elif isinstance(gp, GaussianProcess):
            gp = _fit_single_GP_MAP(gp, n_tries, theta0, method, **kwargs)
        elif LibGPGPU.gpu_usable() and isinstance(gp, GaussianProcessGPU):
            gp = _fit_single_GPGPU_MAP(gp, n_tries, theta0, method, **kwargs)
        elif LibGPGPU.gpu_usable() and isinstance(gp, MultiOutputGP_GPU):
            gp = _fit_MOGPGPU_MAP(gp, n_tries, theta0, method, **kwargs)
        else:
            raise TypeError("single arg to fit_GP_MAP must be a GaussianProcess or MultiOutputGP instance")
    elif len(args) < 2:
        raise TypeError("missing required inputs/targets arrays to GaussianProcess")
    else:
        gp_kwargs = {}
        for key in ["mean", "kernel", "priors", "nugget", "inputdict", "use_patsy"]:
            if key in kwargs:
                gp_kwargs[key] = kwargs[key]
                del kwargs[key]
        try:
            gp = GaussianProcess(*args, **gp_kwargs)
            gp = _fit_single_GP_MAP(gp, n_tries, theta0, method, **kwargs)
        except AssertionError:
            try:
                gp = MultiOutputGP(*args, **gp_kwargs)
                gp = _fit_MOGP_MAP(gp, n_tries, theta0, method, **kwargs)
            except AssertionError:
                raise ValueError("Bad values for *args in fit_GP_MAP")
    if isinstance(gp, GaussianProcessBase):
        if (isinstance(gp, GaussianProcess) and gp.theta.get_data() is None) or \
            (isinstance(gp, GaussianProcessGPU) and gp.theta is None):
            raise RuntimeError("GP fitting failed")
    else:
        if len(gp.get_indices_not_fit()) > 0:
            failure_string = "Fitting failed for emulators {}".format(
                    gp.get_indices_not_fit())
            if skip_failures:
                print(failure_string)
            else:
                raise RuntimeError(failure_string)
    return gp

def _fit_single_GPGPU_MAP(gp, n_tries=15, theta0=None, method='L-BFGS-B', **kwargs):
    """Fit hyperparameters using MAP for a single GP in the C++/CUDA implementation
       The optimization is done in C++, this Python function is a wrapper for that.
        Returns a single GP object that has its hyperparameters fit to the MAP value.
    """
    if method not in ["L-BFGS", "L-BFGS-B"]:
        raise NotImplementedError("Unknown method for optimizer - only L-BFGS implemented for GPU")
    n_tries = int(n_tries)
    assert n_tries > 0, "number of attempts must be positive"
    if theta0 is None or len(theta0)==0:
        theta0=np.array([])
    LibGPGPU.fit_GP_MAP(gp._densegp_gpu, n_tries, theta0)
    return gp

def _fit_MOGPGPU_MAP(gp, n_tries=15, theta0=None, method='L-BFGS-B', **kwargs):
    """Fit hyperparameters using MAP for a multi-output GP in the C++/CUDA implementation
       The optimization is done in C++, this Python function is a wrapper for that.
        Returns an MOGP object that has its hyperparameters fit to the MAP value.
    """
    if method not in ["L-BFGS", "L-BFGS-B"]:
        raise NotImplementedError("Unknown method for optimizer - only L-BFGS implemented for GPU")
    n_tries = int(n_tries)
    assert n_tries > 0, "number of attempts must be positive"
    if theta0 is None or len(theta0)==0:
        theta0=np.array([])
    LibGPGPU.fit_GP_MAP(gp._mogp_gpu, n_tries, theta0)
    return gp

def _fit_single_GP_MAP(gp, n_tries=15, theta0=None, method='L-BFGS-B', **kwargs):
    """Fit hyperparameters using MAP for a single GP

    Returns a single GP object that has its hyperparameters fit to the
    MAP value.  Accepts keyword arguments passed to scipy's
    minimization routine.

    """

    assert isinstance(gp, GaussianProcessBase)
    n_tries = int(n_tries)
    assert n_tries > 0, "number of attempts must be positive"

    np.seterr(divide = 'raise', over = 'raise', invalid = 'raise')

    logpost_values = []
    theta_values = []

    for i in range(n_tries):
        if i == 0 and not theta0 is None:
            theta = np.array(theta0)
            assert theta.shape == (gp.n_params,), "theta0 must be a 1D array with length n_params"
        else:
            theta = gp.priors.sample()
        try:
            min_dict = minimize(gp.logposterior, theta, method = method,
                                jac = gp.logpost_deriv, options = kwargs)

            min_theta = min_dict['x']
            min_logpost = min_dict['fun']

            logpost_values.append(min_logpost)
            theta_values.append(min_theta)
        except LinAlgError:
            print("Matrix not positive definite, skipping this iteration")
        except FloatingPointError:
            print("Floating point error in optimization routine, skipping this iteration")

    if len(logpost_values) == 0:
        print("Minimization routine failed to return a value")
        gp.theta = None
    else:
        logpost_values = np.array(logpost_values)
        idx = np.argmin(logpost_values)

        gp.fit(theta_values[idx])

    return gp

def _fit_single_GP_MAP_bound(gp, theta0, n_tries, method, **kwargs):
    "fitting function accepting theta0 as an argument for parallelization"

    return _fit_single_GP_MAP(gp, n_tries=n_tries, theta0=theta0, method=method, **kwargs)

def _fit_MOGP_MAP(gp, n_tries=15, theta0=None, method='L-BFGS-B',
                  refit=False, **kwargs):
    """Fit hyperparameters using MAP for multiple GPs in parallel

    Uses Python Multiprocessing to fit GPs in parallel by calling the
    above routine for a single GP for each of the emulators in the
    MOGP class. Returns a MultiOutputGP object where all emulators
    have been fit to the MAP value.

    Routine only fits GPs that have not previously been fit (indicated
    by the hyperparameters being set to ``None`` by default. This can
    be overridden by passing ``refit=True``.

    Accepts a ``processes`` argument (integer or None) as a keyword to
    control the number of subprocesses used to fit the individual GPs
    in parallel. Must be positive. Default is ``None``.

    """

    assert isinstance(gp, MultiOutputGP)

    try:
        processes = kwargs['processes']
        del kwargs['processes']
    except KeyError:
        processes = None

    n_tries = int(n_tries)
    assert n_tries > 0, "n_tries must be a positive integer"

    if theta0 is None:
        theta0 = [ None ]*gp.n_emulators
    else:
        if isinstance(theta0, np.ndarray):
            if theta0.ndim == 1:
                theta0 = [theta0]*gp.n_emulators
            else:
                assert theta0.ndim == 2, "theta0 must be a 1D or 2D array"
                assert theta0.shape[0] == gp.n_emulators, "bad shape for fitting starting points"
        elif isinstance(theta0, list):
            assert len(theta0) == gp.n_emulators, "theta0 must be a list of length n_emulators"

    if not processes is None:
        processes = int(processes)
        assert processes > 0, "number of processes must be positive"

    if refit:
        emulators_to_fit = gp.emulators
        indices_to_fit = list(range(len(gp.emulators)))
        thetavals = theta0
    else:
        indices_to_fit = gp.get_indices_not_fit()
        emulators_to_fit = gp.get_emulators_not_fit()
        thetavals = [ theta0[idx] for idx in indices_to_fit]

    if platform.system() == "Windows":

        fit_MOGP = [fit_GP_MAP(emulator, n_tries=n_tries, theta0=t0, method=method, **kwargs)
                    for (emulator, t0) in zip(emulators_to_fit, thetavals)]
    else:
        with Pool(processes) as p:
            fit_MOGP = p.starmap(partial(_fit_single_GP_MAP_bound, n_tries=n_tries, method=method, **kwargs),
                                 [(emulator, t0) for (emulator, t0) in zip(emulators_to_fit, thetavals)])

    for (idx, em) in zip(indices_to_fit, fit_MOGP):
        gp.emulators[idx] = em

    return gp

max_stars_repo_name: EXAUQ/mogp-emulator

max_issues_repo_name: EXAUQ/mogp-emulator

max_forks_repo_name: EXAUQ/mogp-emulator

reason: import numpy,import scipy,from scipy


text:
"""A module containing numpy implementations of the Keras layers
"""
import numpy as np

class Layer(object):
    """Dummy layer base class for numpy layers
    """
    def __init__(self, *args, **kwargs):
        self._is_built = False

    def build(self, inputs):
        pass

    def call(self, inputs):
        raise NotImplementedError

    def __call__(self, *args, **kwargs):
        if not self._is_built:
            self.build(*args, **kwargs)
            self._is_built = True
        return self.call(*args, **kwargs)

#
# Basic math functions
#
class OneHot(Layer):
    """One-hot atomic number layer

    Converts a list of atomic numbers to one-hot vectors

    Input: atomic numbers (batch, atoms)
    Output: one-hot atomic number (batch, atoms, atomic_number)
    """
    def __init__(self,
                 max_atomic_number,
                 **kwargs):
        # Parameters
        self.max_atomic_number = max_atomic_number

        super(OneHot, self).__init__(**kwargs)

    def call(self, inputs):
        atomic_numbers = inputs
        eye_mat = np.eye(self.max_atomic_number + 1)
        return eye_mat[atomic_numbers]

class DistanceMatrix(Layer):
    """
    Distance matrix layer

    Expands Cartesian coordinates into a distance matrix.

    Input: coordinates (..., atoms, 3)
    Output: distance matrix (..., atoms, atoms)
    """
    def call(self, positions):
        # `positions` should be Cartesian coordinates of shape
        #    (..., atoms, 3)
        v1 = np.expand_dims(positions, axis=-2)
        v2 = np.expand_dims(positions, axis=-3)

        sum_squares = np.sum(np.square(v2 - v1), axis=-1)
        return np.sqrt(sum_squares)

class AngleTensor(Layer):
    """
    Angle tensor layer

    Calculate the tensor containing all possible angles between
    cartesian positions (between [0, 180] degrees)

    The returned tensor is indexed as follows:
        the first n_points axis is the central point forming the angle
        the latter two n_points axes are the two other points

    The tensor should be symmetric about a permutation of the last two indices.

    Input: coordinates (..., atoms, 3)
    Output: distance matrix (..., atoms, atoms)
    """
    def __init__(self, deg=False, eps=1e-10):
        super(AngleTensor, self).__init__()
        self.eps = eps
        self.deg = deg

    def call(self, inputs):
        # `positions` should be Cartesian coordinates of shape
        #    (..., atoms, 3)
        positions, magnitudes = inputs

        v1 = np.expand_dims(positions, axis=-2)
        v2 = np.expand_dims(positions, axis=-3)

        diff = v2 - v1
        magnitude_products = (
            np.expand_dims(magnitudes, axis=-1)
            * np.expand_dims(magnitudes, axis=-2)
        )
        dot_prod = np.sum(
            np.expand_dims(diff, axis=-2) * np.expand_dims(diff, axis=-3),
            axis=-1
        )
        # Avoids division by zero
        magnitude_products[magnitude_products < self.eps] = 1.0

        # Calculate the angles
        angles = np.arccos(dot_prod / magnitude_products)

        # Zero invalid values
        n_positions = positions.shape[-2]
        mask1 = np.reshape(1 - np.eye(n_positions, dtype=int), (1, n_positions, n_positions))
        mask2 = np.reshape(1 - np.eye(n_positions, dtype=int), (n_positions, 1, n_positions))
        mask3 = np.reshape(1 - np.eye(n_positions, dtype=int), (n_positions, n_positions, 1))
        mask = mask1 * mask2 * mask3

        angles *= mask

        # Convert to degrees if asked
        if self.deg:
            angles *= 180 / np.pi

        return angles

#
# Kernel functions
#
class KernelBasis(Layer):
    """Expand tensor using kernel of width=width, spacing=spacing,
    starting at min_value ending at max_value (inclusive if endpoint=True).

    Input: tensor (batch, atoms, [atoms, [atoms...])
    Output: tensor expanded into kernel basis set (batch, atoms, [atoms, [atoms...]], n_gaussians)

    Args:
        min_value (float, optional): minimum value
        max_value (float, optional): maximum value (non-inclusive)
        width (float, optional): width of kernel functions
        spacing (float, optional): spacing between kernel functions
        self_thresh (float, optional): value below which a distance is
            considered to be a self interaction (i.e. zero)
        include_self_interactions (bool, optional): whether or not to include
            self-interactions (i.e. distance is zero)
    """
    def __init__(self, min_value=-1, max_value=9, width=0.2,
                 spacing=0.2, self_thresh=1e-5, include_self_interactions=True,
                 endpoint=False, **kwargs):
        super(KernelBasis, self).__init__(**kwargs)
        self._n_centers = int(np.ceil((max_value - min_value) / spacing))
        self.min_value = min_value
        self.max_value = max_value
        self.spacing = spacing
        self.width = width
        self.self_thresh = self_thresh
        self.include_self_interactions = include_self_interactions
        self.endpoint = endpoint
        self.mu = None

    def call(self, in_tensor):
        in_tensor = np.expand_dims(in_tensor, -1)
        self.mu = np.linspace(self.min_value, self.max_value, self._n_centers,
                              endpoint=self.endpoint)

        mu_prefix_shape = tuple([1 for _ in range(len(in_tensor.shape) - 1)])
        mu = np.reshape(self.mu, mu_prefix_shape + (-1,))
        values = self.kernel_func(in_tensor, mu)

        if not self.include_self_interactions:
            mask = (in_tensor >= self.self_thresh).astype(in_tensor.dtype)
            values *= mask

        return values

    def kernel_func(self, inputs, centres):
        raise NotImplementedError

class GaussianBasis(KernelBasis):
    """Expand distance matrix into Gaussians of width=width, spacing=spacing,
    starting at min_value ending at max_value (inclusive if endpoint=True).

        -(x - u)^2
    exp(----------)
        2 * w^2

    where: u is linspace(min_value, max_value, ceil((max_value - min_value) / width))
           w is width

    Input: distance_matrix (batch, atoms, atoms)
    Output: distance_matrix expanded into Gaussian basis set (batch, atoms, atoms, n_centres)

    Args:
        min_value (float, optional): minimum value
        max_value (float, optional): maximum value (non-inclusive)
        width (float, optional): width of Gaussians
        spacing (float, optional): spacing between Gaussians
        self_thresh (float, optional): value below which a distance is
            considered to be a self interaction (i.e. zero)
        include_self_interactions (bool, optional): whether or not to include
            self-interactions (i.e. distance is zero)
                (batch, atoms, atoms, n_gaussians)
    """
    def kernel_func(self, inputs, centres):
        gamma = -0.5 / (self.width  ** 2)
        return np.exp(gamma * np.square(inputs - centres))

class TriangularBasis(KernelBasis):
    """Expand distance matrix into triangles of width=width, spacing=spacing,
    starting at min_value ending at max_value (inclusive if endpoint=True).

          1
    1 - -----|x - u|
         2w

    where: u is linspace(min_value, max_value, ceil((max_value - min_value) / width))
           w is width

    Input: distance_matrix (batch, atoms, atoms)
    Output: distance_matrix expanded into Triangular basis set (batch, atoms, atoms, n_centres)

    Args:
        min_value (float, optional): minimum value
        max_value (float, optional): maximum value (non-inclusive)
        width (float, optional): width of triangles
        spacing (float, optional): spacing between triangles
        self_thresh (float, optional): value below which a distance is
            considered to be a self interaction (i.e. zero)
        include_self_interactions (bool, optional): whether or not to include
            self-interactions (i.e. distance is zero)
                (batch, atoms, atoms, n_triangles)
    """
    def kernel_func(self, inputs, centres):
        gamma = -0.5 / self.width
        return 1 + gamma * np.abs(inputs - centres)

#
# Cutoff functions
#
class CutoffLayer(Layer):
    """Base layer for cutoff functions.

    Applies a cutoff function to the expanded distance matrix

    Inputs:
        distance_matrix (batch, atoms, atoms)
        basis_functions (batch, atoms, atoms, n_centres)
    Output: basis_functions with cutoff function multiplied (batch, atoms, atoms, n_centres)
    """
    def __init__(self, cutoff, **kwargs):
        super(CutoffLayer, self).__init__(**kwargs)
        self.cutoff = cutoff

    def call(self, inputs):
        distance_matrix, basis_functions = inputs

        cutoffs = self.cutoff_function(distance_matrix)
        cutoffs = np.expand_dims(cutoffs, axis=-1)

        return basis_functions * cutoffs

    def cutoff_function(self, distance_matrix):
        """Function responsible for the cutoff. It should also return zeros
        for anything greater than the cutoff.

        Args:
            distance_matrix (Tensor): the distance matrix tensor
        """
        raise NotImplementedError

class CosineCutoff(CutoffLayer):
    """The cosine cutoff originally proposed by Behler et al. for ACSFs.
    """
    def cutoff_function(self, distance_matrix):
        cos_component = 0.5 * (1 + np.cos(np.pi * distance_matrix / self.cutoff))
        return np.where(
            distance_matrix <= self.cutoff,
            cos_component,
            np.zeros_like(distance_matrix)
        )

class TanhCutoff(CutoffLayer):
    """Alternate tanh^3 cutoff function mentioned in some of the ACSF papers.
    """
    def cutoff_function(self, distance_matrix):
        normalization_factor = 1.0 / (np.tanh(1.0) ** 3)
        tanh_component = (np.tanh(1.0 - (distance_matrix / self.cutoff))) ** 3
        return np.where(
            distance_matrix <= self.cutoff,
            normalization_factor * tanh_component,
            np.zeros_like(distance_matrix)
        )

class LongTanhCutoff(CutoffLayer):
    """Custom tanh cutoff function that keeps symmetry functions relatively unscaled
    longer than the previously proposed tanh function
    """
    def cutoff_function(self, distance_matrix):
        normalization_factor = 1.0 / (np.tanh(float(self.cutoff)) ** 3)
        tanh_component = (np.tanh(self.cutoff - distance_matrix)) ** 3
        return np.where(
            distance_matrix <= self.cutoff,
            normalization_factor * tanh_component,
            np.zeros_like(distance_matrix)
        )

#
# Atom-related functions
#
class AtomicNumberBasis(Layer):
    """Expands Gaussian matrix into the one-hot atomic numbers basis

    Inputs:
        one_hot_numbers  (batch, atoms, max_atomic_number + 1)
        gaussians_matrix  (batch, atoms, atoms, n_gaussians)
    Output:
        gaussians_atom_matrix  (batch, atoms, atoms, n_gaussians, max_atomic_number + 1)
    """
    def __init__(self, zero_dummy_atoms=False, **kwargs):
        super(AtomicNumberBasis, self).__init__(**kwargs)
        self.zero_dummy_atoms = zero_dummy_atoms

    def call(self, inputs):
        one_hot_numbers, gaussian_mat = inputs

        gaussian_mat = np.expand_dims(gaussian_mat, axis=-1)
        if self.zero_dummy_atoms:
            mask = np.eye(one_hot_numbers.shape[-1], dtype=gaussian_mat.dtype)
            mask[0] = 0
            one_hot_numbers = np.tensordot(one_hot_numbers, mask, axes=1)
        one_hot_numbers = np.expand_dims(one_hot_numbers, axis=1)
        one_hot_numbers = np.expand_dims(one_hot_numbers, axis=3)
        return gaussian_mat * one_hot_numbers

#
# Normalization-related layers
#
class Unstandardization(Layer):
    """
    Offsets energies by mean and standard deviation (optionally, per-atom)

    `mu` and `sigma` both follow the following:
        If the value is a scalar, apply it equally to all properties
        and all types of atoms

        If the value is a vector, each component corresponds to an
        output property. It is expanded into a matrix where the
        first axis shape is 1. It then follows the matrix rules.

        If the value is a matrix, rows correspond to types of atoms and
        columns correspond to properties.

            If there is only one row, then the row vector applies to every
            type of atom equally.

            If there is one column, then the scalars are applied to every
            property equally.

            If there is a single scalar, then it is treated as a scalar.

    Inputs: the inputs to this layer depend on whether or not mu and sigma
            are given as a single scalar or per atom type.

        If scalar:
            atomic_props  (batch, atoms, energies)
        If per type:
            one_hot_atomic_numbers (batch, atoms, atomic_number)
            atomic_props  (batch, atoms, energies)
    Output: atomic_props  (batch, atoms, energies)

    Attributes:
        mu (float, list, or np.ndarray): the mean values by which
            to offset the inputs to this layer
        sigma (float, list, or np.ndarray): the standard deviation
            values by which to scale the inputs to this layer
    """
    def __init__(self, mu, sigma, trainable=False, per_type=None, use_float64=False, **kwargs):
        super(Unstandardization, self).__init__(trainable=trainable, **kwargs)
        self.init_mu = mu
        self.init_sigma = sigma
        self.use_float64 = use_float64
        self.dtype = 'float64' if use_float64 else 'float32'

        self.mu = np.asanyarray(self.init_mu)
        self.sigma = np.asanyarray(self.init_sigma)

        self.per_type = len(self.mu.shape) > 0 or per_type

    @staticmethod
    def expand_ones_to_shape(arr, shape):
        if len(arr.shape) == 0:
            arr = arr.reshape((1 ,1))
        if 1 in arr.shape:
            tile_shape = tuple(shape[i] if arr.shape[i] == 1 else 1
                               for i in range(len(shape)))
            arr = np.tile(arr, tile_shape)
        if arr.shape != shape:
            raise ValueError('the arrays were not of the right shape: '
                             'expected %s but was %s' % (shape, arr.shape))
        return arr

    def build(self, inputs):
        # If mu is given as a vector, assume it applies to all atoms
        if len(self.mu.shape) == 1:
            self.mu = np.expand_dims(self.mu, axis=0)
        if len(self.sigma.shape) == 1:
            self.sigma = np.expand_dims(self.sigma, axis=0)

        if self.per_type:
            one_hot_atomic_numbers, atomic_props = inputs
            one_hot_atomic_numbers = one_hot_atomic_numbers.shape
            atomic_props = atomic_props.shape

            w_shape = (one_hot_atomic_numbers[-1], atomic_props[-1])

            self.mu = self.expand_ones_to_shape(self.mu, w_shape)
            self.sigma = self.expand_ones_to_shape(self.sigma, w_shape)
        else:
            w_shape = self.mu.shape

        self.mu = self.mu.astype(self.dtype)
        self.sigma = self.sigma.astype(self.dtype)

    def call(self, inputs):
        # `atomic_props` should be of shape (batch, atoms, energies)

        # If mu and sigma are given per atom type, need atomic numbers
        # to know how to apply them. Otherwise, just energies is enough.
        if self.per_type or isinstance(inputs, (list, tuple)):
            one_hot_atomic_numbers, atomic_props = inputs
        else:
            atomic_props = inputs
        atomic_props = atomic_props.astype(self.dtype)
        one_hot_atomic_numbers = one_hot_atomic_numbers.astype(self.dtype)

        if self.per_type:
            atomic_props *= np.tensordot(one_hot_atomic_numbers, self.sigma, axes=1)
            atomic_props += np.tensordot(one_hot_atomic_numbers, self.mu, axes=1)
        else:
            atomic_props *= self.sigma
            atomic_props += self.mu

        return atomic_props

class AtomRefOffset(Unstandardization):
    """
    Offsets energies by per-atom reference properties.
    Simpler case of Unstandardization.

    Inputs: one_hot_atomic_numbers (batch, atoms, atomic_number)
            atomic_props  (batch, atoms, energies)
    Output: atomic_props  (batch, atoms, energies)

    Attributes:
        atom_ref (list or np.ndarray): atom references of
            shape (atomic_number, n_props)
    """
    def __init__(self, atom_ref=None, add_offset=True, **kwargs):
        self.add_offset = add_offset
        self.atom_ref = atom_ref

        if not self.add_offset:
            self.atom_ref = -self.atom_ref

        kwargs.pop('per_type', None)
        super(AtomRefOffset, self).__init__(
            mu=self.atom_ref,
            sigma=kwargs.pop('sigma', 1.0),
            per_type=True,
            **kwargs
        )

#
# Dummy atom-related layers
#
class DummyAtomMasking(Layer):
    """
    Masks dummy atoms (atomic number = 0 by default) with zeros

    Inputs: atomic_numbers
                Either or both in this order:
                    atomic_numbers  (batch, atoms)
                or
                    one_hot_atomic_numbers  (batch, atoms, atomic_number)
            value  (batch, atoms, ...)
    Output: value with zeroes for dummy atoms  (batch, atoms, ...)

    Args:
        atom_axes (int or iterable of int): axes to which to apply
            the masking

    Keyword Args:
        dummy_index (int): the index to mask (default: 0)
        invert_mask (bool): if True, zeroes all but the desired index rather
            than zeroeing the desired index
    """
    def __init__(self, atom_axes=1, **kwargs):
        self.invert_mask = kwargs.pop('invert_mask', False)
        self.dummy_index = kwargs.pop('dummy_index', 0)
        super(DummyAtomMasking, self).__init__(**kwargs)
        if isinstance(atom_axes, int):
            atom_axes = [atom_axes]
        elif isinstance(atom_axes, tuple):
            atom_axes = list(atom_axes)
        self.atom_axes = atom_axes

    def call(self, inputs):
        # `value` should be of shape (batch, atoms, ...)
        one_hot_atomic_numbers, value = inputs
        atomic_numbers = np.argmax(one_hot_atomic_numbers,
                                   axis=-1)

        # Form the mask that removes dummy atoms (atomic number = dummy_index)
        if self.invert_mask:
            selection_mask = atomic_numbers == self.dummy_index
        else:
            selection_mask = atomic_numbers != self.dummy_index
        selection_mask = selection_mask.astype(value.dtype)

        for axis in self.atom_axes:
            mask = selection_mask
            for _ in range(axis - 1):
                mask = np.expand_dims(mask, axis=1)
            # Add one since K.int_shape does not return batch dim
            while len(value.shape) != len(mask.shape):
                mask = np.expand_dims(mask, axis=-1)

            # Zeros the energies of dummy atoms
            value *= mask
        return value

class SelectAtoms(DummyAtomMasking):
    """Selects atoms of a certain type.

    Convenience function around DummyAtomMasking to have it
    zero all other atoms rather than just the dummy index.
    """
    def __init__(self, *args, **kwargs):
        atom_index = kwargs.pop('atom_index', kwargs.pop('dummy_index', 0))
        kwargs['invert_mask'] = True
        super(SelectAtoms, self).__init__(
            dummy_index=atom_index,
            *args, **kwargs
        )

max_stars_repo_name: UPEIChemistry/atomic-images

max_issues_repo_name: UPEIChemistry/atomic-images

max_forks_repo_name: UPEIChemistry/atomic-images


text:
import numpy as np

from PuzzleLib.Backend import gpuarray
from PuzzleLib.Modules.Module import ModuleError, Module

class Concat(Module):
	def __init__(self, axis, name=None):
		super().__init__(name)
		self.registerBlueprint(locals())

		self.axis = axis
		self.sections = None

	def updateData(self, data):
		self.sections = [d.shape[self.axis] for d in data]
		self.data = gpuarray.concatenate(data, axis=self.axis)

	def updateGrad(self, grad):
		self.grad = gpuarray.split(grad, self.sections, axis=self.axis)

	def checkDataShape(self, shapes):
		for i, shape in enumerate(shapes[1:]):
			if not shape[:self.axis] + shape[self.axis + 1:] == shapes[0][:self.axis] + shapes[0][self.axis + 1:]:
				raise ModuleError(
					"Shape %d is inconsistent with initial shape (checking %s, init is %s)" % (i, shape, shapes[0])
				)

	def dataShapeFrom(self, shapes):
		concatDim = 0
		for shape in shapes:
			concatDim += shape[self.axis]

		shape = shapes[0][:self.axis] + (concatDim, ) + shapes[0][self.axis + 1:]
		return shape

	def checkGradShape(self, shape):
		concatDim = 0
		for sec in self.sections:
			concatDim += sec

		gradShape = self.data.shape[:self.axis] + (concatDim, ) + self.data.shape[self.axis+1:]
		if gradShape != shape:
			raise ModuleError("Expected grad shape %s (given %s)" % (gradShape, shape))

	def gradShapeFrom(self, shape):
		shapes = []
		for sec in self.sections:
			shapes.append(shape[:self.axis] + (sec, ) + shape[self.axis + 1:])

		return shapes

	def calcMode(self, T):
		dtypes = {dtype for dtype, _ in gpuarray.dtypesSupported()}

		if T not in dtypes:
			raise ModuleError("Unsupported dtype %s" % T)

		self.calctype = T

def unittest():
	alongBatchAxis()
	alongDataAxis()

def alongBatchAxis():
	data = []
	for _ in range(3):
		data.append(gpuarray.to_gpu(np.random.randn(np.random.randint(low=5, high=10), 10, 5, 3).astype(np.float32)))

	concat = Concat(axis=0)
	concat(data)

	hostOutData = np.concatenate([d.get() for d in data], axis=0)
	assert np.allclose(hostOutData, concat.data.get())

	grad = gpuarray.to_gpu(np.random.randn(*hostOutData.shape).astype(np.float32))
	concat.backward(grad)

	stride = 0
	hostInGrad = []
	for i in range(len(data)):
		hostInGrad.append(grad.get()[stride:stride + data[i].shape[0], :])
		stride += data[i].shape[0]

	assert all([np.allclose(hostInGrad[i], concat.grad[i].get()) for i in range(len(data))])

def alongDataAxis():
	data = []
	for _ in range(3):
		data.append(gpuarray.to_gpu(np.random.randn(10, np.random.randint(low=4, high=8), 4, 5).astype(np.float32)))

	concat = Concat(axis=1)
	concat(data)

	hostOutData = np.concatenate([d.get() for d in data], axis=1)
	assert np.allclose(hostOutData, concat.data.get())

	grad = gpuarray.to_gpu(np.random.randn(*hostOutData.shape).astype(np.float32))
	concat.backward(grad)

	stride = 0
	hostInGrad = []
	for i in range(len(data)):
		hostInGrad.append(grad.get()[:, stride:stride + data[i].shape[1]])
		stride += data[i].shape[1]

	assert all([np.allclose(hostInGrad[i], concat.grad[i].get()) for i in range(len(data))])

if __name__ == "__main__":
	unittest()

max_stars_repo_name: EmilPi/PuzzleLib

max_issues_repo_name: EmilPi/PuzzleLib

max_forks_repo_name: EmilPi/PuzzleLib


text:
#!/usr/bin/env python3
# -*- coding: utf-8 -*-

import pandas as pd
import numpy as np
import json
from scipy import interpolate
import logging
from tespy.networks import load_network
from tespy.tools import logger
from tespy.connections import Ref
from tespy.tools.helpers import TESPyNetworkError
from tespy.tools import document_model

logger.define_logging(
    log_path=True, log_version=True, screen_level=logging.WARNING, file_level=logging.WARNING
)

class PowerPlantModel:
    """
    Class for a power plant model within a coupled loop.

    Parameters
    ----------
    path : str
        Path to exported tespy model
    """

    def __init__(self, path):

        self.path = path
        self.load_tespy_model()

    def load_tespy_model(self):

        self.nw = load_network(self.path)
        self.nw.set_attr(iterinfo=False)
        self.nw.solve('design', init_only=True)

    def get_parameters(self, **kwargs):

        result = kwargs.copy()
        if "Connections" in kwargs:
            for c, params in kwargs["Connections"].items():
                for param in params:
                    result["Connections"][c][param] = self.nw.get_conn(c).get_attr(param)

        if "Components" in kwargs:
            for c, params in kwargs["Components"].items():
                for param in params:
                    result["Components"][c][param] = self.nw.get_comp(c).get_attr(param)

    def get_single_parameter(self, obj, label, parameter):
        if obj == 'Components':
            return self.nw.get_comp(label).get_attr(parameter).val
        elif obj == 'Connections':
            return self.nw.get_conn(label).get_attr(parameter).val

    def set_parameters(self, **kwargs):

        if "Connections" in kwargs:
            for c, params in kwargs["Connections"].items():
                self.nw.get_conn(c).set_attr(**params)

        if "Components" in kwargs:
            for c, params in kwargs["Components"].items():
                self.nw.get_comp(c).set_attr(**params)

    def set_single_parameter(self, obj, label, parameter, value):
        if obj == 'Components':
            self.nw.get_comp(label).set_attr({parameter: value})
        elif obj == 'Connections':
            self.nw.get_conn(label).get_attr({parameter: value})

    def run(self, **kwargs):

        self.set_parameters(**kwargs)

        self.solved = False
        try:
            self.nw.solve("offdesign", design_path=self.path)
            if self.nw.res[-1] >= 1e-3 or self.nw.lin_dep:
                self.nw.solve(
                    "offdesign", init_only=True, init_path=self.path,
                    design_path=self.path
                )
            else:
                if any(self.nw.results['HeatExchanger']['Q'] > 0):
                    self.solved = False
                else:
                    self.solved = True
        except ValueError as e:
            print(e)
            print("The iteration failed with the following exception: " + e)
            self.nw.lin_dep = True
            self.nw.solve(
                "offdesign", init_only=True, init_path=self.path,
                design_path=self.path
            )

reason: import numpy,from scipy


text:
import numpy as np

turn = np.zeros(14, dtype=np.float32)
turnTmp = "00000000000000"
ans = 4
if ans < 14:
    turnTmp = turnTmp[:ans] + "1"+turnTmp[ans+1:]
print(turnTmp)

for i in range(14):
    turn[i] = float(turnTmp[i])
print(turn)

max_stars_repo_name: km-t/dcpython

max_issues_repo_name: km-t/dcpython

max_forks_repo_name: km-t/dcpython


text:
"""
Please note, this code is only for python 3+. If you are using python 2+, please modify the code accordingly.
"""
import tensorflow as tf
import numpy as np

# create data
x_data = np.random.rand(100).astype(np.float32)
y_data = x_data*0.1 + 0.3

### create tensorflow structure start ###

### create tensorflow structure end ###
          # Very important

for step in range(201):
    pass


text:
import time
import math
import dolfin as df
import numpy as np

"""
This script works with dolfin 1.0.1.

In dolfin 1.1.0, the ordering of degrees of freedom does not correspond
to geometrical structures such as vertices anymore. This has to be
investigated further. For now, the re-ordering can be turned off using
the commented line below.

"""
# df.parameters.reorder_dofs_serial = False # for dolfin 1.1.0

"""
The goal is to compute the values of f(r, t) = sin(x) * sin(t) on a mesh for
different points in time. The cost of updating f a couple of times is measured.

"""

L = math.pi / 2; n = 10;
mesh = df.Box(0, 0, 0, L, L, L, n, n, n)
#mesh = df.BoxMesh(0, 0, 0, L, L, L, n, n, n) # dor dolfin 1.1.0
ts = np.linspace(0, math.pi / 2, 100)

# dolfin expression code

S = df.FunctionSpace(mesh, "CG", 1)
expr = df.Expression("sin(x[0]) * sin(t)", t = 0.0)

start = time.time()
for t in ts:
    expr.t = t
    f_dolfin = df.interpolate(expr, S)
t_dolfin = time.time() - start
print "Time needed for dolfin expression: {:.2g} s.".format(t_dolfin)

# explicit loop code

f_loop = np.empty(mesh.num_vertices())

start = time.time()
for t in ts:
    for i, (x, y, z) in enumerate(mesh.coordinates()):
        f_loop[i] = math.sin(x) * math.sin(t)
t_loop = time.time() - start
print "Time needed for loop: {:.2g} s.".format(t_loop)

# compute the ratio

ratio = t_dolfin / t_loop
if ratio >= 1:
    print "Looping over numpy array is {:.2g} times faster than interpolating dolfin expression.".format(ratio)
else:
    print "Interpolating the dolfin expression is {:.2g} times faster than looping over numpy array.".format(ratio)

# check correctness

if False:
    print "\nMesh coordinates:\n", mesh.coordinates()
    print "Function was sin(x) * sin(t)."
    print "Dolfin vector at t=pi/2:\n\t", f_dolfin.vector().array()
    print "Numpy vector at t=pi/2:\n\t", f_loop

max_diff = np.max(np.abs(f_dolfin.vector().array() - f_loop))
assert max_diff < 1e-14, "Maximum difference is d = {:.2g}.".format(max_diff)

max_stars_repo_name: davidcortesortuno/finmag

max_issues_repo_name: davidcortesortuno/finmag

max_forks_repo_name: davidcortesortuno/finmag


text:
"""This is a script to test MetaReplayBuffer."""
import numpy as np

from metarl.replay_buffer.meta_replay_buffer import MetaReplayBuffer
from metarl.sampler import PEARLSampler
from metarl.tf.envs import TfEnv
from metarl.torch.policies import DeterministicMLPPolicy
from tests.fixtures.envs.dummy import DummyBoxEnv

def test_all():
    """Test all methods."""
    env = TfEnv(DummyBoxEnv())
    policy = DeterministicMLPPolicy(env_spec=env.spec, hidden_sizes=(32, 32))
    max_path_length = 10
    max_samples = 50
    max_trajs = 50

    sampler = PEARLSampler(env, policy, max_path_length)

    paths, _ = sampler.obtain_samples(max_samples=max_samples,
                                      max_trajs=max_trajs,
                                      accum_context=False)

    replay_buffer = MetaReplayBuffer(100, env.observation_space.low.size,
                                     env.action_space.low.size)

    i = 0
    for path in paths:
        replay_buffer.add_path(path)
        i += max_path_length
        assert replay_buffer.size() == i

    replay_buffer.clear()
    assert replay_buffer.size() == 0

    for path in paths:
        replay_buffer.add_path(path)

    batch_size = 3
    indices = np.random.randint(0, replay_buffer.size(), batch_size)
    out = replay_buffer.sample_data(indices)
    assert len(out['observations']) == batch_size
    assert len(out['actions']) == batch_size
    assert len(out['rewards']) == batch_size
    assert len(out['terminals']) == batch_size
    assert len(out['next_observations']) == batch_size

    batch_size = 10
    out = replay_buffer.random_batch(batch_size)
    assert len(out['observations']) == batch_size
    assert len(out['actions']) == batch_size
    assert len(out['rewards']) == batch_size
    assert len(out['terminals']) == batch_size
    assert len(out['next_observations']) == batch_size

    out = replay_buffer.random_sequence(batch_size)
    assert len(out['observations']) == batch_size
    assert len(out['actions']) == batch_size
    assert len(out['rewards']) == batch_size
    assert len(out['terminals']) == batch_size
    assert len(out['next_observations']) == batch_size


text:
# Copyright (c) 2021 PaddlePaddle Authors. All Rights Reserved.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
"""
temp script
"""
import sys
sys.path.append("..")
import paddle
import numpy as np
from adversary import Adversary
from attacks.gradient_method import FGSM
from attacks.cw import CW_L2
from attacks.saliency import JSMA
from models.whitebox import PaddleWhiteBoxModel

from classifier.definednet import transform_eval, TowerNet
model_0 = TowerNet(3, 10, wide_scale=1)
model_1 = TowerNet(3, 10, wide_scale=2)

# set fgsm attack configuration
fgsm_attack_config = {"norm_ord": np.inf, "epsilons": 0.003, "epsilon_steps": 1, "steps": 1}
paddle_model = PaddleWhiteBoxModel(
    [model_0, model_1],   # ensemble two models
    [1, 1.8],             # dictate weight
    paddle.nn.CrossEntropyLoss(),
    (-3, 3),
    channel_axis=3,
    nb_classes=10)

# FGSM attack, init attack with the ensembled model
# attack = FGSM(paddle_model)
attack = CW_L2(paddle_model, learning_rate=0.01)
# attack = JSMA(paddle_model)

cifar10_test = paddle.vision.datasets.Cifar10(mode='test', transform=transform_eval)
test_loader = paddle.io.DataLoader(cifar10_test, batch_size=1)

data = test_loader().next()
img = data[0][0]
label = data[1]

# init adversary status
adversary = Adversary(img.numpy(), int(label))
target = np.random.randint(paddle_model.num_classes())
while label == target:
    target = np.random.randint(paddle_model.num_classes())
print(label, target)
adversary.set_status(is_targeted_attack=True, target_label=target)

# launch attack
# adversary = attack(adversary, **fgsm_attack_config)
adversary = attack(adversary, attack_iterations=100, verbose=True)

if adversary.is_successful():
    original_img = adversary.original
    adversarial_img = adversary.adversarial_example
    print("Attack succeeded.")
else:
    print("Attack failed.")

max_stars_repo_name: csdongxian/PaddleSleeve

max_issues_repo_name: csdongxian/PaddleSleeve

max_forks_repo_name: csdongxian/PaddleSleeve


text:
from pyns.phy import PHYLayer
import scipy.special
import math

class BPSK(PHYLayer):
    def __init__(self, threshold, bandwidth = 12500):
        super().__init__(threshold, bandwidth, 1)

    def compute_ber(self, ebn0, is_log=False):
        # use in db
        if is_log: ebn0 = 10 ** (ebn0 / 10)
        return 0.5 * scipy.special.erfc(math.sqrt(ebn0))

class QPSK(PHYLayer):
    def __init__(self, threshold, bandwidth = 12500):
        super().__init__(threshold, bandwidth, 2)

    def compute_ber(self, ebn0, is_log=False):
        # use in db
        if is_log: ebn0 = 10 ** (ebn0 / 10)
        # QPSK union bound
        return 2 * scipy.special.erfc(math.sqrt(2 * ebn0)) + \
                scipy.special.erfc(math.sqrt(4 * ebn0))

max_stars_repo_name: Kuree/pyns

max_issues_repo_name: Kuree/pyns

max_forks_repo_name: Kuree/pyns

reason: import scipy


text:
import numpy as np

from .vector import shift

class synchronization():
    def __init__(self, param):
        self.halfpreamble = param.halfpreamble
        self.full_pream_len = len(self.halfpreamble) * 2
        self.Ncp = param.Ncp
        self.Ncs = param.Ncs
        self.K = param.K
        self.M = param.M
        self.N = param.K * param.M

    def detect_preamble_starts(self, data):
        self.half_peam_len = len(self.halfpreamble)
        metric = self._calc_metric(data)

        peak_locs = self._find_metric_peaks(metric)
        #peak_locs = self._find_single_peak(metric)
        preamble_starts = peak_locs - int(self.half_peam_len / 2)
        # print(preamble_starts)
        return preamble_starts

    def _calc_metric(self, data):
        cross_corr = np.correlate(data, self.halfpreamble, mode="same")
        metric = cross_corr + shift(cross_corr, -self.half_peam_len)
        metric = (metric.real)**2 + (metric.imag)**2
        autocorr = self._calc_moving_autocorr(data)
        return metric * autocorr

    def _calc_threshold(self, metric):
        half_peam_len = len(self.halfpreamble)
        threshold = np.empty(len(metric))
        metric = np.pad(metric, (half_peam_len, 0), 'constant')
        for i in range(half_peam_len + 1, len(metric)):
            window = metric[i - half_peam_len:i]
            window[np.where(window > threshold[i - half_peam_len - 1])] /= 2
            threshold[i - half_peam_len] = np.sum(window)
        threshold = shift(threshold / 2, half_peam_len //
                              4, mode="same", fill_value=threshold[-1] / 2)
        threshold[np.where(threshold < max(metric) / 100)] = max(metric) / 100
        return threshold

    def _find_metric_peaks(self, metric):
        threshold = self._calc_threshold(metric)
        # step 1: find peaks
        locs = np.where(
            (metric > shift(metric, 1)) & (metric > shift(metric, -1, fill_value=metric[-1])) &
            (metric > threshold)
        )[0]
        # step 2: find max peak in window with length of fullpreamble
        last_peak = metric[locs[0]]
        last_loc = locs[0]
        locs_out = np.array([], int)
        for l in locs:
            if ((l - last_loc) > (self.half_peam_len + 5)):
                locs_out = np.append(locs_out, last_loc)
                last_peak = 0.0
            if last_peak < metric[l]:
                last_peak = metric[l]
                last_loc = l
        locs_out = np.append(locs_out, last_loc)
        return locs_out

    def _find_single_peak(self, metric):
        threshold = self._calc_threshold(metric)
        loc = np.argmax(metric)
        if metric[loc] > threshold[loc]:
            return np.array([loc])
        else:
            return np.array([])

    def _calc_moving_autocorr(self, data):
        half_peam_len = len(self.halfpreamble)
        Ncp_cs = self.Ncp + self.Ncs
        autocorr = data * shift(data, half_peam_len)
        autocorr_metric = np.empty_like(autocorr)
        for i in range(half_peam_len, len(autocorr)):
            autocorr_metric[i -
                            half_peam_len] = np.sum(autocorr[i - half_peam_len:i])
        autocorr_metric = np.abs(autocorr_metric)
        autocorr_out = np.empty_like(autocorr_metric)
        for i in range(Ncp_cs - 2, len(autocorr)):
            autocorr_out[i - Ncp_cs -
                         2] = np.sum(autocorr_metric[i - Ncp_cs - 2:i])
        # shift(autocorr_out,Ncp_cs+half_peam_len)
        return autocorr_out / max(autocorr_out)

max_stars_repo_name: jdemel/GFDM-PHY-Reference

max_issues_repo_name: jdemel/GFDM-PHY-Reference

max_forks_repo_name: jdemel/GFDM-PHY-Reference


text:
#!/usr/bin/env python
#
# Copyright (C) 2009 Rosen Diankov
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
from __future__ import with_statement # for python 2.5
PKG = 'posedetectiondb' # this package name
NAME = 'GatherDetectionResults'

import roslib; roslib.load_manifest(PKG)
import os, sys, time, string, threading
from optparse import OptionParser
import numpy # nice to be able to explicitly call some functions
from numpy import *

import rospy
import posedetection_msgs.msg

from openravepy import *

class VisibilityModel(metaclass.AutoReloader):
    def __init__(self,measurements=None,filename=None,symmetricplane=None,kinbodyfile=None):
        if measurements is not None:
            self.rawmeasurements = measurements

        if symmetricplane is not None:
            dists = dot(symmetricplane[0:3],transpose(self.rawmeasurements))+symmetricplane[3]
            self.rawmeasurements = r_[self.rawmeasurements,self.rawmeasurements-dot(reshape(2.0*dists,(dists.shape[0],1)),reshape(symmetricplane[0:3],(1,3)))]

        self.trimesh = None
        if kinbodyfile is not None:
            self.env = Environment()
            self.orobj = self.env.ReadKinBodyXMLFile(kinbodyfile)
            if self.orobj is None:
                raise ValueError('failed to open %s openrave file'%kinbodyfile)
            self.env.AddKinBody(self.orobj)
            self.trimesh = self.env.Triangulate(self.orobj)
        self.measurements = self.rawmeasurements
    def CreateReducedModel(self,bandwidth=0.04,bandthresh=0.01,neighthresh=0.01,showdata=False,savefile=None):
        self.measurements,indices = self.Prune(self.rawmeasurements,100,neighthresh**2,1)
        uniformpoints,dists,pointscale = self.UniformlySampleSpace(bandwidth,bandthresh)
        if showdata:
            from enthought.mayavi import mlab
            mlab.figure(1,fgcolor=(0,0,0), bgcolor=(1,1,1))
            src = mlab.pipeline.scalar_field(dists)
            mlab.pipeline.iso_surface(src,contours=[0.01],opacity=0.1)
            mlab.pipeline.volume(mlab.pipeline.scalar_field(dists*500))
            v = pointscale[0]*self.trimesh.vertices+pointscale[1]
            mlab.triangular_mesh(v[:,0],v[:,1],v[:,2],self.trimesh.indices,color=(0,0,0.5))
        if savefile is None:
            savefile = self.getfilename()
        print('saving measurements to %s'%savefile)
        mkdir_recursive(os.path.split(savefile)[0])
        savetxt(savefile,uniformpoints,'%f')
    def getfilename(self):
        return os.path.join(self.env.GetHomeDirectory(),'kinbody.'+self.orobj.GetKinematicsGeometryHash(),'visibility.txt')

    def UniformlySampleSpace(self,bandwidth,bandthresh,delta=0.02):
        maxradius = sqrt(max(sum(self.measurements**2,1)))
        nsteps = floor(maxradius/delta)
        X,Y,Z = mgrid[-nsteps:nsteps,-nsteps:nsteps,-nsteps:nsteps]
        allpoints = c_[X.flat,Y.flat,Z.flat]*delta/bandwidth
        sampleinds = flatnonzero(sum(allpoints**2,1)<(maxradius/bandwidth)**2)
        samplepoints = allpoints[sampleinds,:]
        kdtree = pyANN.KDTree(self.measurements/bandwidth)
        sampledists = zeros(samplepoints.shape[0])
        goodpoints = []
        for i in xrange(samplepoints.shape[0]):
            neighs,dists,kball = kdtree.kFRSearchArray(samplepoints[i:(i+1),:],5.0**2,32,0.0001)
            sampledists[i] = sum(exp(-dists[neighs>=0]))
        uniformpoints = samplepoints[sampledists>bandthresh,:]*bandwidth
        alldists = zeros(prod(X.shape))
        alldists[sampleinds] = sampledists
        return uniformpoints,reshape(alldists,X.shape),array((1.0/delta,nsteps))

    def Prune(self,rawposes, nsize, thresh2, neighsize,giveupiters=100):
        """rawposes is Nx7"""
        iter = 1
        poses = array(rawposes)
        indices = range(poses.shape[0])
        N = poses.shape[0]
        nochange=0
        while N > nsize:
            ind = numpy.random.randint(N)
            g = poses[ind,:]
            # check g's neighbors
            d = sum((poses[0:N,:] - tile(g, (N,1)))**2,1)
            neigh = sum(d < thresh2)
            if neigh > neighsize:
                # move to the last pose and resize
                poses[ind,:] = poses[N-1,:]
                indices[ind] = indices[N-1]
                nochange=0
                N -= 1
            nochange += 1
            iter += 1
            if iter > 5000 or nochange > giveupiters:
                break
        return poses[0:N,:],indices[0:N]

class OpenRAVEVisualizer(metaclass.AutoReloader):
    def __init__(self,kinbodyfile,automaticadd=True,measurementsfilename=None):
        self.automaticadd = automaticadd
        self.orenv = Environment()
        self.orenv.SetViewer('qtcoin')
        self.orobj = self.orenv.ReadKinBodyXMLFile(kinbodyfile)
        if self.orobj is None:
            raise ValueError('failed to open %s openrave file'%kinbodyfile)
        self.orenv.AddKinBody(self.orobj)
        self.objab = self.orobj.ComputeAABB()
        self.Tcamera = None
        self.lck = threading.Lock()
        self.camerahandle = None
        self.measurementhandles = []
        self.measurements = []
        if measurementsfilename is not None:
            f = open(measurementsfilename,'r')
            for l in f.readlines():
                m = array([string.atof(s) for s in string.split(l)])
                self.drawmeasurement(m)
                self.measurements.append(m)
            f.close()

        self.sub_objdet = rospy.Subscriber("ObjectDetection", posedetection_msgs.msg.ObjectDetection,self.objdetcb, queue_size=1)
        rospy.init_node(NAME, anonymous=True)#,disable_signals=False)

    def __del__(self):
        self.sub_objdet.unregister()
        self.orenv.Destroy()

    def objdetcb(self,msg):
        newcamerahandle = None
        if len(msg.objects) > 0:
            q = msg.objects[0].pose.orientation
            t = msg.objects[0].pose.position
            with self.lck:
                self.Tcamera = linalg.inv(matrixFromPose([q.w,q.x,q.y,q.z,t.x,t.y,t.z]))

            # draw in openrave environment
            sx = 0.02
            sy = 0.02
            camlocalpts = array(((sx,sy,0.05),(-sx,sy,0.05),(-sx,-sy,0.05),(sx,-sy,0.05),(sx,sy,0.05),(0,0,0),(-sx,sy,0.05),(0,0,0),(-sx,-sy,0.05),(0,0,0),(sx,-sy,0.05)))
            campts = dot(self.Tcamera,r_[transpose(camlocalpts),ones((1,camlocalpts.shape[0]))])
            self.camerahandle = self.orenv.drawlinestrip(transpose(campts[0:3,:]),2,array([0,0,1]))
        if self.automaticadd:
            self.addmeasurement()

    def addmeasurement(self):
        self.lck.acquire()
        Tcamera = self.Tcamera
        self.lck.release()
        dist = dot(Tcamera[0:3,2],self.objab.pos()-Tcamera[0:3,3:4])
        m = Tcamera[0:3,2]*dist
        self.drawmeasurement(m)
        self.measurements.append(m)
        print('num measurements %d'%len(self.measurements))

    def drawmeasurement(self,m):
        dist = sqrt(sum(m**2))
        dir = m/dist
        p = self.objab.pos()-dist*dir
        self.measurementhandles.append(self.orenv.drawlinestrip(transpose(c_[p-0.02*dir,p+0.02*dir]),4,array((1,min(1,dist/1.0),0))))

    def savemeasurements(self,filename):
        self.lck.acquire()
        print('saving measurements to %s'%filename)
        savetxt(filename,self.measurements,'%f')
        self.lck.release()

if __name__=='__main__':
    parser = OptionParser(description='Gather object detection transformations and filter and display them.')
    parser.add_option('--kinbodyfile',
                      action="store",type='string',dest='kinbodyfile',
                      help='OpenRAVE object file that represents the incoming object pose. Updates are show in the openrave window')
    parser.add_option('-s','--single',
                      action="store_true",dest='single',default=False,
                      help='If set, will wait for user input in order to add a measurement')
    parser.add_option('-f','--savefile',
                      action="store",dest="filename",
                      help='If specified, will save all recorded measurements to this file at exit time')
    parser.add_option('-m','--measurements',
                      action="store",dest="measurements",default=None,
                      help='If specified, will start with the current measurements file')
    (options, args) = parser.parse_args()
    if not options.kinbodyfile:
        print('Error: Need to specify an openrave kinbody file')
        sys.exit(1)

    visualizer = OpenRAVEVisualizer(options.kinbodyfile,measurementsfilename=options.measurements,automaticadd=not options.single)
    while True:
        cmd = raw_input('Enter command (q-quit and save,c-capture): ');
        if cmd == 'q':
            break
        elif cmd == 'c' and options.single:
            print('adding measurement')
            visualizer.addmeasurement()
        else:
            print('bad command',cmd)
    if options.filename:
        visualizer.savemeasurements(options.filename)

def test():
    "rosrun posedetectiondb GatherDetectionResults.py --kinbodyfile=scenes/cereal_frootloops.kinbody.xml -f test.txt ObjectDetection:=/CerealDetection"
    import GatherDetectionResults
    self = GatherDetectionResults.VisibilityModel(measurements=loadtxt('uroncha_raw.txt'),symmetricplane=array([1.0,0,0,0]),kinbodyfile='scenes/uroncha.kinbody.xml')
    self = GatherDetectionResults.VisibilityModel(measurements=loadtxt('peachjuice_raw.txt'),kinbodyfile='scenes/peachjuice.kinbody.xml')
    self.CreateReducedModel(bandwidth=0.03,bandthresh=0.02)

reason: import numpy,from numpy


text:
"""Test module for progress and reaction rate methods."""

import os
from flamespeed.test import env
from flamespeed import chemkin
import numpy as np

# =================
# TEST PARAM
# =================

def test_threebody_param():
    """Test progress rate values."""
    a = chemkin.ReactionRate()
    a.read_XML('./data/rxn_ThreeBody.xml', verify_integrity=True, convert_units=True)
    test = [[{'Type': 'Kooij', 'A': 38000000000.0, 'b': -2.0, 'E': 0.0, 'name': None},
            {'Type': 'efficiencies', 'H2': '2.5', 'H2O': '12.0', 'default': '1.0'}],
            [{'Type': 'Kooij', 'A': 48000000000.0, 'b': -5.0, 'E': 4.1840000000000002, 'name': None},
            {'Type': 'efficiencies', 'H2': '2.5', 'H2O': '12.0', 'default': '1.0'}
            ]]
    assert (a.k_params_nonelementary == test)

def test_threebodytroe_param():
    """Test progress rate values."""
    a = chemkin.ReactionRate()
    a.read_XML('./data/rxn_TroeFalloffThreeBody.xml', verify_integrity=True, convert_units=True)
    test = [[{'Type': 'Kooij', 'A': 636600000.0, 'b': -1.72, 'E': 2195.7631999999999, 'name': 'k0'},
                {'Type': 'Kooij', 'A': 1475000.0, 'b': 0.6, 'E': 0.0, 'name': None},
                {'Type': 'Troe', 'alpha': 0.8, 't1': 1e+30, 't2': 1e+30, 't3': 1e+30},
                {'Type': 'efficiencies', 'H2': '2.0', 'H2O': '11.0', 'O2': '0.78', 'default': '1.0'}
                ]]
    assert (a.k_params_nonelementary == test)

def test_duplicate_param():
    """Test progress rate values."""
    a = chemkin.ReactionRate()
    a.read_XML('./data/rxn_Duplicate.xml', verify_integrity=True, convert_units=True)
    test = [[{'Type': 'Arrhenius', 'A': 420000000.0, 'E': 50132.688000000002},
                {'Type': 'Arrhenius', 'A': 130000.0, 'E': -6816.9912000000004}]]
    assert (a.k_params_nonelementary == test)

max_stars_repo_name: bkleyn/flamespeed

max_issues_repo_name: bkleyn/flamespeed

max_forks_repo_name: bkleyn/flamespeed


text:
import torch
import numpy as np
from torch.utils.data import Dataset

from dataset import DataTaker

class InsuranceAnswerDataset(Dataset):

    def __init__(self, data_type='train', dataset_size=100, negative_size=0, full_flag=False, full_size=0):
        """
        :param data_type: 'train', 'valid' or 'test'
        :param dataset_size: relate to the dataset size, but not equivalent
        :param negative_size: negative answer number for one positive answer
        """
        self.data_list = []
        self.dataset_size = dataset_size
        self.negative_size = negative_size
        if full_flag:
            self.dataset_size = 100
        dt = DataTaker(dataset_size=self.dataset_size, full_flag=full_flag, full_size=full_size)
        if data_type == 'train':
            self.q_list = dt.read_train()
        elif data_type == 'test':
            self.q_list = dt.read_test()
        elif data_type == 'valid':
            self.q_list = dt.read_valid()
        self.a_list = dt.read_answer()
        self.max_q_length = dt.question_max_len
        self.max_a_length = dt.answer_max_len
        key_id = 0
        for line in self.q_list:
            for gt_id in line['ground_truth']:
                temp_data = []
                question = np.pad(line['question'], (0, self.max_q_length - len(line['question'])), 'constant',
                                  constant_values=0)
                temp_data.append([question, len(line['question']), key_id])
                answer = np.pad(self.a_list[gt_id], (0, self.max_a_length - len(self.a_list[gt_id])), 'constant',
                                constant_values=0)
                temp_data.append([torch.LongTensor(answer), len(self.a_list[gt_id]), key_id])
                nt_id_list = line['negative_pool']
                np.random.shuffle(nt_id_list)
                nt_list = []
                for nt_id in nt_id_list[:self.negative_size]:
                    answer = np.pad(self.a_list[nt_id], (0, self.max_a_length - len(self.a_list[nt_id])), 'constant',
                                    constant_values=0)
                    nt_list.append([torch.LongTensor(answer), len(self.a_list[nt_id]), key_id])
                key_id += 1
                temp_data.append(nt_list)
                del nt_list
                self.data_list.append(temp_data)
                del temp_data

    def __getitem__(self, item):
        data = self.data_list[item]
        return data, 1

    def __len__(self):
        return len(self.data_list)

max_stars_repo_name: zs-liu/Pytorch-AS

max_issues_repo_name: zs-liu/Pytorch-AS

max_forks_repo_name: zs-liu/Pytorch-AS