text:
from __future__ import division
import PIL
# from matplotlib.colors import ListedColormap, LinearSegmentedColormap
from matplotlib import cm
import matplotlib.pyplot as plt
import seaborn as sn

# import shutil
import numpy as np
import torch
# from pathlib import Path
# import datetime
# from collections import OrderedDict
import matplotlib as mpl
mpl.use('Agg')  # No x-server

device = torch.device(
    "cuda") if torch.cuda.is_available() else torch.device("cpu")

# Seaborn settings
# sns.set(style=\"whitegrid\", rc={\"font.size\":8,\"axes.titlesize\":8,\"axes.labelsize\":5})
sn.set(style="whitegrid", font_scale=1.5)
sn.set_palette("bright", n_colors=4, color_codes=True)

# def high_res_colormap(low_res_cmap, resolution=1000, max_value=1):
#     # Construct the list colormap, with interpolated values for higer resolution
#     # For a linear segmented colormap, you can just specify the number of point in
#     # cm.get_cmap(name, lutsize) with the parameter lutsize
#     x = np.linspace(0, 1, low_res_cmap.N)
#     low_res = low_res_cmap(x)
#     new_x = np.linspace(0, max_value, resolution)
#     high_res = np.stack([np.interp(new_x, x, low_res[:, i])
#                          for i in range(low_res.shape[1])], axis=1)
#     return ListedColormap(high_res)

# def opencv_rainbow(resolution=1000):
#     # Construct the opencv equivalent of Rainbow
#     opencv_rainbow_data = (
#         (0.000, (1.00, 0.00, 0.00)),
#         (0.400, (1.00, 1.00, 0.00)),
#         (0.600, (0.00, 1.00, 0.00)),
#         (0.800, (0.00, 0.00, 1.00)),
#         (1.000, (0.60, 0.00, 1.00))
#     )

#     return LinearSegmentedColormap.from_list('opencv_rainbow', opencv_rainbow_data, resolution)

# COLORMAPS = {'rainbow': opencv_rainbow(),
#              'magma': high_res_colormap(cm.get_cmap('magma')),
#              'bone': cm.get_cmap('bone', 1000)}

imagenet_mean = torch.Tensor([0.485, 0.456, 0.406])
imagenet_std = torch.Tensor([0.229, 0.224, 0.225])

# imagenet_mean_grey = torch.Tensor([0.445])
# imagenet_std_grey = torch.Tensor([0.269])

def tensor2array(tensor, max_value=None, colormap='magma'):
    #tensor = tensor.detach().cpu()
    if max_value is None:
        max_value = tensor.abs().max()  # .item()
    if tensor.ndimension() == 2 or tensor.size(0) == 1:
        # if is_grey:  # grey image
        #     array = imagenet_mean_grey[:, None, None].to(tensor.device) + \
        #         tensor*imagenet_std_grey[:, None, None].to(tensor.device)
        #     array = array.detach().cpu().numpy()
        # else:  # Depth image
        # norm_array = tensor.squeeze().numpy()/max_value
        norm_array = tensor.squeeze()/max_value
        norm_array = norm_array.detach().cpu().numpy()
        # array = COLORMAPS[colormap](norm_array).astype(np.float32)
        array = cm.get_cmap(colormap)(norm_array).astype(np.float32)
        array = array.transpose(2, 0, 1)

    elif tensor.ndimension() == 3:
        assert(tensor.size(0) == 3)
        # array = 119.4501 + tensor.numpy()*6.5258
        # array = imagenet_mean[:, None, None].to(tensor.device) + \
        #     tensor*imagenet_std[:, None, None].to(tensor.device)
        array = tensor.detach().cpu().numpy()
    return array

def traj2Img(pred_xyz):
    """Make an image from the 2D plot of a given trajectory.

    Args:
        pred_xyz (torch.Tensor): Trajectory to plot. Shape: [N,3]

    Returns:
        PIL.Image: Image of the trajectory plot
    """

    pred_xyz = pred_xyz.detach().cpu()
    pred_x = pred_xyz[:, 0].numpy()
    pred_y = pred_xyz[:, 1].numpy()

    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.plot(pred_x, pred_y)
    fig.canvas.draw()

    img_plt = PIL.Image.frombytes(
        'RGB', fig.canvas.get_width_height(), fig.canvas.tostring_rgb())
    img_plt = np.asarray(img_plt)  # CHW
    img_plt = img_plt.transpose(2, 0, 1)  # HWC
    return img_plt

def traj2Fig(pred_xyz, axes=[0, 1]):
    """Make `matplotlib.pyplot.figure` from the 2D plot of a given trajectory.

    Args:
        pred_xyz (torch.Tensor): Trajectory to plot. Shape: [N,3]

    Returns:
        matplotlib.pyplot.figure: Figure of the trajectory plot
    """

    pred_xyz = pred_xyz.cpu()

    fig = plt.figure()
    ax = plt.gca()
    ax.plot(pred_xyz[:, axes[0]], pred_xyz[:, axes[1]])
    # fig.canvas.draw()

    return fig

def traj2Fig_withgt(pred_xyz, gt_xyz, axes=[0, 1]):
    """Make `matplotlib.pyplot.figure` from the 2D plot of a given trajectory.

    Args:
        pred_xyz (torch.Tensor): Trajectory to plot. Shape: [N,3]
        gt_xyz (torch.Tensor): Trajectory to plot. Shape: [N,3]

    Returns:
        matplotlib.pyplot.figure: Figure of the trajectory plot
    """

    pred_xyz = pred_xyz.cpu()
    gt_xyz = gt_xyz.cpu()

    fig = plt.figure()
    ax = plt.gca()
    ax.plot(pred_xyz[:, axes[0]], pred_xyz[:, axes[0]], label='Prediction')
    ax.plot(gt_xyz[:, 0], gt_xyz[:, 1], label='Ground-truth')
    ax.legend()
    # fig.canvas.draw()

    return fig

def save_traj_plots(results_dir, f_pred_xyz, b_pred_xyz, axes=[0, 1], plt_prefix=''):
    """Plot and save trajectories to a file.

    Args:
        results_dir (str): Directory to save the plots.
        f_pred_xyz (torch.Tensor): Predicted xyz values of a trajectory.
        b_pred_xyz (torch.Tensor): Predicted inverse xyz values of a trajectory.
        axes (list, optional): Which axes to use in 2D plot. Defaults to [0, 1].
    """
    plt.figure(figsize=(16, 8))
    plt.subplot(1, 2, 1)
    ax = sn.lineplot(x=f_pred_xyz[:, axes[0]].cpu().numpy(
    ), y=f_pred_xyz[:, axes[1]].cpu().numpy(), sort=False)
    ax.set(title='Forward Trajectory', xlabel='X (m)', ylabel='Y (m)')
    plt.subplot(1, 2, 2)
    ax = sn.lineplot(x=b_pred_xyz[:, axes[0]].cpu().numpy(
    ), y=b_pred_xyz[:, axes[1]].cpu().numpy(), sort=False, markers=True)
    ax.set(title='Backward Trajectory', xlabel='X (m)', ylabel='Y (m)')

    # Save fig
    plt.tight_layout()
    plt.savefig(results_dir/(plt_prefix+'pred_nogt.pdf'),
                bbox_inches='tight', pad_inches=0)
    plt.savefig(results_dir/(plt_prefix+'pred_nogt.png'),
                bbox_inches='tight', pad_inches=0)

def save_traj_plots_with_gt(results_dir, pred_xyz, gt, plt_prefix=''):
    """Plot and save trajectories together with ground-truth to a file.

    Args:
        results_dir (str): Directory to save the plots.
        pred_xyz (torch.Tensor): Predicted xyz values of a trajectory.
        gt (torch.Tensor): Ground truth predictions as [4,4] transformation matrices.
    """
    gt_xyz = gt[:, :3, 3].cpu().numpy()
    # np_pred = 0.5*gt_xyz + 0.5*pred_xyz[0].cpu().numpy()
    np_pred = pred_xyz[0].cpu()
    fig, ax = plt.subplots(figsize=(8, 8))
    sn.lineplot(x=np_pred[:, 0], y=np_pred[:, 1],
                sort=False, ax=ax, label='Ours')
    sn.lineplot(x=gt_xyz[:, 0], y=gt_xyz[:, 1],
                sort=False, ax=ax, label='Ground-truth')

    ax.set(xlabel='X (m)', ylabel='Y (m)')

    # Save fig
    plt.tight_layout()
    plt.savefig(results_dir/(plt_prefix+'pred_with_gt.pdf'),
                bbox_inches='tight', pad_inches=0)
    plt.savefig(results_dir/(plt_prefix+'pred_with_gt.png'),
                bbox_inches='tight', pad_inches=0)

# def save_checkpoint_optim(save_path, optim_state, prefix='optim', epoch='', filename='checkpoint.pth.tar'):
#     torch.save(optim_state, save_path /
#                '{}_{}{}'.format(prefix, epoch, filename))

# def save_checkpoint(save_path, masknet_state, posenet_state, epoch='', filename='checkpoint.pth.tar'):
#     file_prefixes = ['radar_posenet']
#     states = [posenet_state]
#     if masknet_state is not None:
#         file_prefixes.append('radar_masknet')
#         states.append(masknet_state)
#     for (prefix, state) in zip(file_prefixes, states):
#         torch.save(state, save_path/'{}_{}{}'.format(prefix, epoch, filename))

#     # if is_best:
#     #     for prefix in file_prefixes:
#     #         shutil.copyfile(save_path/'{}_{}'.format(prefix, filename),
#     #                         save_path/'{}_best.pth.tar'.format(prefix))

# def save_checkpoint_mono(save_path, dispnet_state, vo_pose_state, fusenet_state, epoch='', filename='checkpoint.pth.tar'):
#     file_prefixes = ['mono_dispnet', 'mono_posenet', 'mono_fusenet']
#     states = [dispnet_state, vo_pose_state, fusenet_state]
#     for (prefix, state) in zip(file_prefixes, states):
#         torch.save(state, save_path/'{}_{}{}'.format(prefix, epoch, filename))

#     # if is_best:
#     #     for prefix in file_prefixes:
#     #         shutil.copyfile(save_path/'{}_{}'.format(prefix, filename),
#     #                         save_path/'{}_model_best.pth.tar'.format(prefix))

def save_checkpoint_list(save_path, state_dicts, file_prefixes, epoch='', filename='checkpoint.pth.tar'):
    """Saves a given list of model dictionaries as checpoints.

    Args:
        save_path (str, pathlib.Path): Directory to save the checkpoints.
        state_dicts (list): List of dictionaries to save. It can contain dictionaries and weights of multiple networks and optimizers.
        file_prefixes (list): Names of the models to be saved as prefixes.
        epoch (str, optional): If given, appended to the filename. Defaults to ''.
        filename (str, optional): Default filename postfix. Defaults to 'checkpoint.pth.tar'.
    """
    for (prefix, state) in zip(file_prefixes, state_dicts):
        torch.save(state, save_path/'{}_{}{}'.format(prefix, epoch, filename))

max_stars_repo_name: yasinalm/gramme

max_issues_repo_name: yasinalm/gramme

max_forks_repo_name: yasinalm/gramme


text:
# -*- coding: utf-8 -*-
"""The Dimension object: attributes and methods."""
import warnings
from copy import deepcopy

import numpy as np

from .labeled import LabeledDimension  # lgtm [py/import-own-module]
from .linear import LinearDimension  # lgtm [py/import-own-module]
from .monotonic import MonotonicDimension  # lgtm [py/import-own-module]
from csdmpy.units import string_to_quantity  # lgtm [py/import-own-module]
from csdmpy.utils import _get_dictionary  # lgtm [py/import-own-module]
from csdmpy.utils import validate  # lgtm [py/import-own-module]

__author__ = "Deepansh J. Srivastava"
__email__ = "srivastava.89@osu.edu"
__all__ = ["Dimension"]

functional_dimension = ["linear"]

DEFAULT_DIM = {
    "type": None,  # valid for all dimension subtypes
    "description": "",  # valid for all dimension subtypes
    "count": None,  # valid for linear subtype
    "increment": None,  # valid for linear subtype
    "labels": None,  # valid for labeled subtype
    "coordinates": None,  # valid for monotonic subtype
    "coordinates_offset": None,  # valid for linear subtype
    "origin_offset": None,  # valid for linear subtype
    "complex_fft": False,  # valid for linear subtype
    "period": None,  # valid for monotonic and linear subtypes
    "quantity_name": None,  # valid for monotonic and linear subtypes
    "label": "",  # valid for all dimension subtypes
    "application": {},  # valid for all dimension subtypes
    "reciprocal": {  # valid for monotonic and linear subtypes
        "increment": None,  # valid for monotonic and linear subtypes
        "coordinates_offset": None,  # valid for monotonic and linear subtypes
        "origin_offset": None,  # valid for monotonic and linear subtypes
        "period": None,  # valid for monotonic and linear subtypes
        "quantity_name": None,  # valid for monotonic and linear subtypes
        "label": "",  # valid for monotonic and linear subtypes
        "description": "",  # valid for monotonic and linear subtypes
        "application": {},  # valid for monotonic and linear subtypes
    },
}

class Dimension:
    """Dimension class.

    An instance of this class describes a dimension of a multi-dimensional system.
    In version 1.0 of the CSD model, there are three subtypes of the Dimension class:

    - :ref:`linearDimension_uml`,
    - :ref:`monotonicDimension_uml`, and
    - :ref:`labeledDimension_uml`.

    **Creating an instance of a dimension object**

    There are two ways of creating a new instance of a Dimension class.

    *From a python dictionary containing valid keywords.*

    .. doctest::

        >>> from csdmpy import Dimension
        >>> dimension_dictionary = {
        ...     'type': 'linear',
        ...     'description': 'test',
        ...     'increment': '5 G',
        ...     'count': 10,
        ...     'coordinates_offset': '10 mT',
        ...     'origin_offset': '10 T'
        ... }
        >>> x = Dimension(dimension_dictionary)

    Here, `dimension_dictionary` is the python dictionary.

    *From valid keyword arguments.*

    .. doctest::

        >>> x = Dimension(type = 'linear',
        ...               description = 'test',
        ...               increment = '5 G',
        ...               count = 10,
        ...               coordinates_offset = '10 mT',
        ...               origin_offset = '10 T')
    """

    __slots__ = ("subtype",)

    def __init__(self, *args, **kwargs):
        """Initialize an instance of Dimension object."""
        default = deepcopy(DEFAULT_DIM)
        default_keys = default.keys()
        input_dict = _get_dictionary(*args, **kwargs)
        input_keys = input_dict.keys()

        if "type" not in input_keys:
            raise KeyError("Missing a required 'type' key from the Dimension object.")

        if "reciprocal" in input_keys:
            input_subkeys = input_dict["reciprocal"].keys()

        _ = [
            [default[key].update({subkey: val[subkey]}) for subkey in input_subkeys]
            if key == "reciprocal"
            else default.update({key: val})
            for key, val in input_dict.items()
            if key in default_keys
        ]
        self.__validate_key_value__(default)

        if default["type"] == "labeled":
            self.subtype = LabeledDimension(**default)

        if default["type"] == "monotonic":
            self.subtype = MonotonicDimension(values=default["coordinates"], **default)

        if default["type"] == "linear":
            self.subtype = self._linear(default)

    @staticmethod
    def __validate_key_value__(default):
        _valid_types = ["monotonic", "linear", "labeled"]

        type_ = default["type"]
        message = (
            f"The value, '{type_}', is invalid for the `type` attribute of the "
            "Dimension object. The allowed values are 'monotonic', 'linear' and "
            "'labeled'."
        )

        if default["type"] not in _valid_types:
            raise ValueError(message)

        if default["type"] == "labeled" and default["labels"] is None:
            raise KeyError(
                "Missing a required `labels` key from the LabeledDimension object."
            )

        if default["type"] == "monotonic" and default["coordinates"] is None:
            raise KeyError(
                "Missing a required `coordinates` key from the MonotonicDimension "
                "object."
            )

    def _linear(self, default):
        """Create and assign a linear dimension."""
        missing_key = ["increment", "count"]

        lst = [item for item in missing_key if default[item] is None]
        if lst != []:
            raise KeyError(
                f"Missing a required `{lst[0]}` key from the LinearDimension object."
            )
        validate(default["count"], "count", int)
        return LinearDimension(**default)

    def __repr__(self):
        """String representation of object."""
        return self.subtype.__repr__()

    def __str__(self):
        """String representation of object."""
        return self.subtype.__str__()

    def __eq__(self, other):
        """Overrides the default implementation."""
        other = other.subtype if isinstance(other, Dimension) else other
        return True if self.subtype == other else False

    def __mul__(self, other):
        """Multiply the Dimension object by a scalar."""
        return self.subtype.__mul__(other)

    def __imul__(self, other):
        """Multiply the Dimension object by a scalar, in-place."""
        return self.subtype.__imul__(other)

    def __truediv__(self, other):
        """Divide the Dimension object by a scalar."""
        return self.subtype.__truediv__(other)

    def __itruediv__(self, other):
        """Divide the Dimension object by a scalar, in-place."""
        return self.subtype.__itruediv__(other)

    # ======================================================================= #
    #                          Dimension Attributes                           #
    # ======================================================================= #
    @property
    def absolute_coordinates(self):
        r"""Absolute coordinates, :math:`\bf X_k^{\rm{abs}}`, along the dimension.

        This attribute is only *valid* for quantitative dimensions, that is,
        `linear` and `monotonic` dimensions. The absolute coordinates are given as

        .. math::

            \mathbf{X}_k^\mathrm{abs} = \mathbf{X}_k + o_k \mathbf{1}

        where :math:`\mathbf{X}_k` are the coordinates along the dimension and
        :math:`o_k` is the :attr:`~csdmpy.Dimension.origin_offset`.
        For example, consider

        .. doctest::

            >>> print(x.origin_offset)
            10.0 T
            >>> print(x.coordinates[:5])
            [100. 105. 110. 115. 120.] G

        then the absolute coordinates are

        .. doctest::

            >>> print(x.absolute_coordinates[:5])
            [100100. 100105. 100110. 100115. 100120.] G

        For `linear` dimensions, the order of the `absolute_coordinates`
        further depend on the value of the
        :attr:`~csdmpy.Dimension.complex_fft` attributes. For
        examples, when the value of the `complex_fft` attribute is True,
        the absolute coordinates are

        .. doctest::

            >>> x.complex_fft = True
            >>> print(x.absolute_coordinates[:5])
            [100075. 100080. 100085. 100090. 100095.] G

        .. testsetup::

            >>> x.complex_fft = False

        Returns:
            A Quantity array of absolute coordinates for quantitative dimensions, `i.e`
            `linear` and `monotonic`.

        Raises:
            AttributeError: For labeled dimensions.
        """
        return self.subtype.absolute_coordinates

    @property
    def application(self):
        """Application metadata dictionary of the dimension object.

        .. doctest::

            >>> print(x.application)
            {}

        The application attribute is where an application can place its own
        metadata as a python dictionary object containing application specific
        metadata, using a reverse domain name notation string as the attribute
        key, for example,

        .. doctest::

            >>> x.application = {
            ...     "com.example.myApp" : {
            ...         "myApp_key": "myApp_metadata"
            ...      }
            ... }
            >>> print(x.application)
            {'com.example.myApp': {'myApp_key': 'myApp_metadata'}}

        Returns:
            A python dictionary containing dimension application metadata.
        """
        return self.subtype.application

    @application.setter
    def application(self, value):
        self.subtype.application = value

    @property
    def axis_label(self):
        r"""Formatted string for displaying label along the dimension axis.

        This attribute is not a part of the original core scientific dataset
        model, however, it is a convenient supplementary attribute that provides
        a formated string ready for labeling dimension axes.
        For quantitative dimensions, this attributes returns a string,
        `label / unit`,  if the `label` is a non-empty string, otherwise,
        `quantity_name / unit`. Here
        :attr:`~csdmpy.Dimension.quantity_name` and
        :attr:`~csdmpy.Dimension.label` are the attributes of the
        :ref:`dim_api` instances, and `unit` is the unit associated with the
        coordinates along the dimension. For examples,

        .. doctest::

            >>> x.label
            'field strength'
            >>> x.axis_label
            'field strength / (G)'

        For `labeled` dimensions, this attribute returns `label`.

        Returns:
            A formated string of label.

        Raises:
            AttributeError: When assigned a value.
        """
        return self.subtype.axis_label

    @property
    def coordinates(self):
        r"""Coordinates, :math:`{\bf X}_k`, along the dimension.

        Example:
            >>> print(x.coordinates)
            [100. 105. 110. 115. 120. 125. 130. 135. 140. 145.] G

        For `linear` dimensions, the order of the `coordinates` also depend on the
        value of the :attr:`~csdmpy.Dimension.complex_fft` attributes.
        For examples, when the value of the `complex_fft` attribute is True,
        the coordinates are

        .. doctest::

            >>> x.complex_fft = True
            >>> print(x.coordinates)
            [ 75.  80.  85.  90.  95. 100. 105. 110. 115. 120.] G

        .. testsetup::

            >>> x.complex_fft = False

        Returns:
            A Quantity array of coordinates for quantitative dimensions, `i.e.` `linear`
            and `monotonic`.

        Returns:
            A Numpy array for labeled dimensions.

        Raises:
            AttributeError: For dimensions with subtype `linear`.
        """
        return self.subtype.coordinates

    @coordinates.setter
    def coordinates(self, value):
        self.subtype.coordinates = value

    @property
    def coords(self):
        """Alias for the `coordinates` attribute."""
        return self.coordinates

    @coords.setter
    def coords(self, value):
        self.coordinates = value

    @property
    def data_structure(self):
        """JSON serialized string describing the Dimension class instance.

        This supplementary attribute is useful for a quick preview of the dimension
        object. The attribute cannot be modified.

        .. doctest::

            >>> print(x.data_structure)
            {
              "type": "linear",
              "count": 10,
              "increment": "5.0 G",
              "coordinates_offset": "10.0 mT",
              "origin_offset": "10.0 T",
              "quantity_name": "magnetic flux density",
              "label": "field strength",
              "description": "This is a test",
              "reciprocal": {
                "quantity_name": "electrical mobility"
              }
            }

        Returns:
            A json serialized string of the dimension object.
        Raises:
            AttributeError: When modified.
        """
        return self.subtype.data_structure

    @property
    def description(self):
        """Brief description of the dimension object.

        The default value is an empty string, ''. The attribute may be
        modified, for example,

        .. doctest::

            >>> print(x.description)
            This is a test

            >>> x.description = 'This is a test dimension.'

        Returns:
            A string of UTF-8 allows characters describing the dimension.

        Raises:
            TypeError: When the assigned value is not a string.
        """
        return self.subtype.description

    @description.setter
    def description(self, value):
        self.subtype.description = value

    @property
    def complex_fft(self):
        """If true, the coordinates are the ordered as the output of a complex fft.

        This attribute is only `valid` for the Dimension instances with `linear`
        subtype.
        The value of this attribute is a boolean specifying if the coordinates along
        the dimension are evaluated as the output of a complex fast Fourier transform
        (FFT) routine.
        For example, consider the following Dimension object,

        .. doctest::

            >>> test = Dimension(type='linear', increment = '1', count = 10)
            >>> test.complex_fft
            False
            >>> print(test.coordinates)
            [0. 1. 2. 3. 4. 5. 6. 7. 8. 9.]

            >>> test.complex_fft = True
            >>> print(test.coordinates)
            [-5. -4. -3. -2. -1.  0.  1.  2.  3.  4.]

        Returns:
            A Boolean.

        Raises:
            TypeError: When the assigned value is not a boolean.
        """
        return self.subtype.complex_fft

    @complex_fft.setter
    def complex_fft(self, value):
        self.subtype.complex_fft = value

    @property
    def increment(self):
        """Increment along a `linear` dimension.

        The attribute is only `valid` for Dimension instances with the subtype
        `linear`. When assigning a value, the dimensionality of the value must
        be consistent with the dimensionality of other members specifying the
        dimension.

        Example:
            >>> print(x.increment)
            5.0 G
            >>> x.increment = "0.1 G"
            >>> print(x.coordinates)
            [100.  100.1 100.2 100.3 100.4 100.5 100.6 100.7 100.8 100.9] G

        Returns:
            A Quantity instance with the increment along the dimension.

        Raises:
            AttributeError: For dimension with subtypes other than `linear`.
            TypeError: When the assigned value is not a string containing a quantity
                       or a Quantity object.
        """
        # .. note:: The sampling interval along a grid dimension and the
        #     respective reciprocal grid dimension follow the Nyquist–Shannon
        #     sampling theorem. Therefore, updating the ``increment``
        #     will automatically trigger an update on its reciprocal
        #     counterpart.
        return self.subtype.increment

    @increment.setter
    def increment(self, value):
        self.subtype.increment = value

    @property
    def coordinates_offset(self):
        r"""Offset corresponding to the zero of the indexes array, :math:`\mathbf{J}_k`.

        When assigning a value, the dimensionality of the value must be consistent with
        the dimensionality of the other members specifying the dimension.

        Example:
            >>> print(x.coordinates_offset)
            10.0 mT
            >>> x.coordinates_offset = "0 T"
            >>> print(x.coordinates)
            [ 0.  5. 10. 15. 20. 25. 30. 35. 40. 45.] G

        The attribute is `invalid` for `labeled` dimensions.

        Returns:
            A Quantity instance with the coordinates offset.

        Raises:
            AttributeError: For `labeled` dimensions.
            TypeError: When the assigned value is not a string containing a quantity
                       or a Quantity object.
        """
        return self.subtype.coordinates_offset

    @coordinates_offset.setter
    def coordinates_offset(self, value):
        self.subtype.coordinates_offset = value

    @property
    def label(self):
        """Label associated with the dimension.

        Example:
            >>> print(x.label)
            field strength
            >>> x.label = 'magnetic field strength'

        Returns:
            A string containing the label.
        Raises:
            TypeError: When the assigned value is not a string.
        """
        return self.subtype.label

    @label.setter
    def label(self, label=""):
        self.subtype.label = label

    @property
    def count(self):
        r"""Number of coordinates, :math:`N_k \ge 1`, along the dimension.

        Example:
            >>> print(x.count)
            10
            >>> x.count = 5

        Returns:
            An Integer specifying the number of coordinates along the dimension.

        Raises:
            TypeError: When the assigned value is not an integer.
        """
        return self.subtype._count

    @count.setter
    def count(self, value):
        self.subtype.count = value

    @property
    def origin_offset(self):
        """Origin offset, :math:`o_k`, along the dimension.

        When assigning a value, the dimensionality of the value must be consistent
        with the dimensionality of other members specifying the dimension.

        Example:
            >>> print(x.origin_offset)
            10.0 T
            >>> x.origin_offset = "1e5 G"

        The origin offset only affect the absolute_coordinates along the dimension.
        This attribute is `invalid` for `labeled` dimensions.

        Returns:
            A Quantity instance with the origin offset.

        Raises:
            AttributeError: For `labeled` dimensions.
            TypeError: When the assigned value is not a string containing a quantity
                       or a Quantity object.
        """
        return self.subtype.origin_offset

    @origin_offset.setter
    def origin_offset(self, value):
        self.subtype.origin_offset = value

    @property
    def period(self):
        """Period of the dimension.

        The default value of the period is infinity, i.e., the dimension is
        non-periodic.

        Example:
            >>> print(x.period)
            inf G
            >>> x.period = '1 T'

        To assign a dimension as non-periodic, one of the following may be
        used,

        .. doctest::

            >>> x.period = '1/0 T'
            >>> x.period = 'infinity µT'
            >>> x.period = '∞ G'

        .. Attention::
            The physical quantity of the period must be consistent with other
            physical quantities specifying the dimension.

        Returns:
            A Quantity instance with the period of the dimension.

        Raises:
            AttributeError: For `labeled` dimensions.
            TypeError: When the assigned value is not a string containing a quantity
                       or a Quantity object.
        """
        return self.subtype.period

    @period.setter
    def period(self, value=None):
        self.subtype.period = value

    @property
    def quantity_name(self):
        """Quantity name associated with the physical quantities specifying dimension.

        The attribute is `invalid` for the labeled dimension.

        .. doctest::

            >>> print(x.quantity_name)
            magnetic flux density

        Returns:
            A string with the `quantity name`.

        Raises:
            AttributeError: For `labeled` dimensions.
            NotImplementedError: When assigning a value.
        """
        return self.subtype.quantity_name

    @quantity_name.setter
    def quantity_name(self, value):
        self.subtype.quantity_name = value

    @property
    def type(self):
        """The dimension subtype.

        There are three *valid* subtypes of Dimension class. The valid
        literals are given by the :ref:`dimObjectSubtype_uml` enumeration.

        .. doctest::

            >>> print(x.type)
            linear

        Returns:
            A string with a valid dimension subtype.

        Raises:
            AttributeError: When the attribute is modified.
        """
        return self.subtype.type

    @property
    def labels(self):
        """Ordered list of labels along the `Labeled` dimension.

        Consider the following labeled dimension,

        .. doctest::

            >>> x2 = Dimension(
            ...         type='labeled',
            ...         labels=['Cu', 'Ag', 'Au']
            ...      )

        then the lables along the labeled dimension are

        .. doctest::

            >>> print(x2.labels)
            ['Cu' 'Ag' 'Au']

        .. note::
            For Labeled dimension, the :attr:`~csdmpy.Dimension.coordinates`
            attribute is an alias of :attr:`~csdmpy.Dimension.labels`
            attribute. For example,

            .. doctest::

                >>> np.all(x2.coordinates == x2.labels)
                True

        In the above example, ``x2`` is an instance of the :ref:`dim_api` class with
        `labeled` subtype.

        Returns:
             A Numpy array with labels along the dimension.

        Raises:
            AttributeError: For dimensions with subtype other than `labeled`.
        """
        return self.coordinates

    @labels.setter
    def labels(self, array):
        self.subtype.labels = array

    @property
    def reciprocal(self):
        """An instance of the ReciprocalDimension class.

        The attributes of ReciprocalDimension class are:
            - coordinates_offset
            - origin_offset
            - period
            - quantity_name
            - label
        where the definition of each attribute is the same as the corresponding
        attribute from the Dimension instance.
        """
        return self.subtype.reciprocal

    # ======================================================================= #
    #                           Dimension Methods                             #
    # ======================================================================= #

    def _copy_metadata(self, obj, copy=False):
        """Copy Dimension metadata"""
        self.subtype._type = obj.subtype._type
        self.subtype._copy_metadata(obj)

    def to_dict(self):
        """Alias to the `dict()` method of the class."""
        return self.dict()

    def dict(self):
        """Return Dimension object as a python dictionary.

        Example:
            >>> x.dict() # doctest: +SKIP
            {'type': 'linear', 'description': 'This is a test', 'count': 10,
            'increment': '5.0 G', 'coordinates_offset': '10.0 mT',
            'origin_offset': '10.0 T', 'quantity_name': 'magnetic flux density',
            'label': 'field strength'}
        """
        return self.subtype.dict()

    def is_quantitative(self):
        """Return True if the dependent variable is quantitative.

        Example:
            >>> x.is_quantitative()
            True
        """
        return self.subtype.is_quantitative()

    def to(self, unit="", equivalencies=None):
        r"""Convert the coordinates along the dimension to the unit, `unit`.

        This method is a wrapper of the `to` method from the
        `Quantity <http://docs.astropy.org/en/stable/api/\
        astropy.units.Quantity.html#astropy.units.Quantity.to>`_ class
        and is only `valid` for physical dimensions.

        Example:
            >>> print(x.coordinates)
            [100. 105. 110. 115. 120. 125. 130. 135. 140. 145.] G
            >>> x.to('mT')
            >>> print(x.coordinates)
            [10.  10.5 11.  11.5 12.  12.5 13.  13.5 14.  14.5] mT

        Args:
            `unit` : A string containing a unit with the same dimensionality as the
                     coordinates along the dimension.

        Raises:
            AttributeError: For `labeled` dimensions.
        """
        self.subtype.to(unit, equivalencies)

    def copy(self):
        """Return a copy of the Dimension object."""
        return deepcopy(self)

    def reciprocal_coordinates(self):
        """Return reciprocal coordinates assuming Nyquist-shannan theorem."""
        return self.subtype.reciprocal_coordinates()

    def reciprocal_increment(self):
        """Return reciprocal increment assuming Nyquist-shannan theorem."""
        return self.subtype.reciprocal_coordinates()

def as_dimension(array, unit="", type=None, **kwargs):
    """Generate and return a Dimension object from a 1D numpy array.

    Args:
        array: A 1D numpy array.
        unit: The unit of the coordinates along the dimension.
        type: The dimension type. Valid values are linear, monotonic, labeled, or
                None. If the value is None, let us decide. The default value is None.
        kwargs: Additional keyword arguments from the Dimension class.

    Example:
        >>> array = np.arange(15)*0.5
        >>> dim_object = cp.as_dimension(array)
        >>> print(dim_object)
        LinearDimension([0.  0.5 1.  1.5 2.  2.5 3.  3.5 4.  4.5 5.  5.5 6.  6.5 7. ])

        >>> array = ['The', 'great', 'circle']
        >>> dim_object = cp.as_dimension(array, label='in the sky')
        >>> print(dim_object)
        LabeledDimension(['The' 'great' 'circle'])
    """
    array = __check_array_for_dimension__(array, type)

    if type is None:
        return _generic_dimensions(array, unit, **kwargs)

    if type == "linear":
        obj = _linear_dimension(array, unit, **kwargs)
        if obj is not None:
            return obj
        raise ValueError("Invalid array for LinearDimension object.")

    if type == "monotonic":
        obj = _monotonic_dimension(array, unit, **kwargs)
        if obj is not None:
            return obj
        raise ValueError("Invalid array for MonotonicDimension object.")

    if type == "labeled":
        if unit != "":
            warnings.warn("Ignoring unit argument for LabeledDimension object.")
        return LabeledDimension(labels=array.tolist(), **kwargs)

def __check_array_for_dimension__(array, type):
    options = [None, "linear", "monotonic", "labeled"]
    if type not in options:
        raise ValueError(f"Invalid value for `type`. Allowed values are {options}.")

    if not isinstance(array, (list, np.ndarray)):
        name = array.__class__.__name__
        raise ValueError(f"Cannot convert {name} to a Dimension object.")

    array = np.asarray(array)
    n = array.ndim
    if n == 1:
        return array

    raise ValueError(f"Cannot convert a {n} dimensional array to a Dimension object.")

def _generic_dimensions(array, unit, className="Dimension", **kwargs):
    """Return a dimension object based on the array coordinates."""
    # labeled
    if str(array.dtype)[:2] in [">U", "<U"]:
        if unit != "":
            warnings.warn("Ignoring unit argument for LabeledDimension.")
        return LabeledDimension(labels=array.tolist(), **kwargs)

    # linear
    obj = _linear_dimension(array, unit, className, **kwargs)
    if obj is not None:
        return obj

    # monotonic
    obj = _monotonic_dimension(array, unit, **kwargs)
    if obj is not None:
        return obj

    raise ValueError("Invalid array for Dimension object.")

def _linear_dimension(array, unit, className="LinearDimension", **kwargs):
    """Return a LinearDimension is array is linear, else None."""
    increment = array[1] - array[0]
    if increment == 0:
        raise ValueError(f"Invalid array for {className} object.")

    if np.allclose(np.diff(array, 1), increment):
        unit = f"({unit})" if str(unit) != "" else ""
        return LinearDimension(
            count=array.size,
            increment=f"{increment} {unit}".strip(),
            coordinates_offset=f"{array[0]} {unit}".strip(),
            **kwargs,
        )

def _monotonic_dimension(array, unit, **kwargs):
    """Return a MonotonicDimension is array is monotonic, else None."""
    if np.all(np.diff(array, 1) > 0) or np.all(np.diff(array, 1) < 0):
        return MonotonicDimension(
            coordinates=array * string_to_quantity(unit), **kwargs
        )

ext: py

max_stars_repo_name: DeepanshS/csdmpy

max_issues_repo_name: deepanshs/csdmpy

max_forks_repo_name: deepanshs/csdmpy

reason: import numpy


text:
import time

from mpmath import *

from numerical_routines import compute_double_pendulum_step_rk4, compute_double_pendulum_step_with_adaptive_step_size_method, get_point_position_mp, get_total_energy_of_pendulum_mp, SimulationAlgorithm, ADAPTIVE_TIME_STEP_METHODS

# Pendulum parameters
gravity = mpf(1)
point1Mass = mpf(1)
point2Mass = mpf(1)
pendulum1Length = mpf(1)
pendulum2Length = mpf(1)
point1AngularVelocity = mpf(0)
point2AngularVelocity = mpf(0)

# Simulation parameters
# mp.prec = 24
mp.prec = 53
# mp.prec = 113

point1Angle = mpf((-3.371910665006095 - -3.396454357612266) / 2 + -3.396454357612266) # center of stable area
point2Angle = mpf((1.925992646191392 - 1.901448953585222) / 2 + 1.901448953585222) # center of stable area
# point1Angle = mpf((-3.39507991082632 - -3.3828080645232346) / 2 + -3.3828080645232346) # edge of stable area
# point2Angle = mpf((1.907830313662826 - 1.9201021599659112) / 2 + 1.9201021599659112) # edge of stable area
# point1Angle = mpf(-3.264360)
# point2Angle = mpf(1.446579)

origin = [mpf(0), mpf(0)]
timeStep = mpf(.01/2**2)
printCurrentStateEveryNSeconds = 2**10
maxTimeToSeeIfPendulumFlips = mpf(2**6)
errorTolerance = mpf(9e-15)
# simulationAlgorithm = SimulationAlgorithm.RK4
# simulationAlgorithm = SimulationAlgorithm.RKF45
# simulationAlgorithm = SimulationAlgorithm.CASH_KARP
# simulationAlgorithm = SimulationAlgorithm.DORMAND_PRINCE_54
# simulationAlgorithm = SimulationAlgorithm.VERNER_65
simulationAlgorithm = SimulationAlgorithm.FEHLBERG_87

# Print simulation info.
print('algorithm = ' + str(simulationAlgorithm.name))
print('binary precision = ' + str(mp.prec))

if simulationAlgorithm in ADAPTIVE_TIME_STEP_METHODS:
    print('error tolerance = ' + str(errorTolerance))

elif simulationAlgorithm is SimulationAlgorithm.RK_4:
    print('time step = ' + str(timeStep))

# Initialization
initialTotalEnergy = get_total_energy_of_pendulum_mp(origin,
                                                     point1Angle, point2Angle,
                                                     point1AngularVelocity, point2AngularVelocity,
                                                     pendulum1Length, pendulum2Length,
                                                     point1Mass, point2Mass,
                                                     gravity)

point1OriginalPosition = get_point_position_mp(origin, point1Angle, pendulum1Length)
elapsedTime = mpf(0)

# Simulate the pendulum.
start = time.time()
timeSinceLastReport = time.time()
while elapsedTime < maxTimeToSeeIfPendulumFlips:
    # RK4
    if simulationAlgorithm is SimulationAlgorithm.RK_4:
        point1Angle, \
        point2Angle, \
        point1AngularVelocity, \
        point2AngularVelocity = compute_double_pendulum_step_rk4(point1Mass, point2Mass,
                                                                 gravity,
                                                                 pendulum1Length, pendulum2Length,
                                                                 point1Angle, point2Angle,
                                                                 point1AngularVelocity, point2AngularVelocity,
                                                                 timeStep)
        elapsedTime += timeStep

    # RKF45
    elif simulationAlgorithm in ADAPTIVE_TIME_STEP_METHODS:
        point1Angle, \
        point2Angle, \
        point1AngularVelocity, \
        point2AngularVelocity, \
        timeStepUsedInCalculation, \
        newTimeStep = compute_double_pendulum_step_with_adaptive_step_size_method(point1Mass, point2Mass,
                                                                                gravity,
                                                                                pendulum1Length, pendulum2Length,
                                                                                point1Angle, point2Angle,
                                                                                point1AngularVelocity, point2AngularVelocity,
                                                                                timeStep,
                                                                                errorTolerance,
                                                                                simulationAlgorithm)

        elapsedTime += timeStepUsedInCalculation
        timeStep = newTimeStep

    # Stop if the pendulum flipped.
    point1CurrentPosition = get_point_position_mp(origin, point1Angle, pendulum1Length)
    if point1CurrentPosition[0]*point1OriginalPosition[0] < 0 and point1CurrentPosition[1] > 0:
        print('the pendulum flipped')
        break
    point1OriginalPosition = point1CurrentPosition

    # Print the current pendulum state at regular intervals.
    if time.time() - timeSinceLastReport > printCurrentStateEveryNSeconds:
        print('simulation time elapsed seconds = ' + str(elapsedTime))
        print('wall time elapsed seconds = ' + str(time.time() - start))
        print('point1Angle = ' + str(point1Angle))
        print('point2Angle = ' + str(point2Angle))
        print('point1AngularVelocity = ' + str(point1AngularVelocity))
        print('point2AngularVelocity = ' + str(point2AngularVelocity))
        print('initial total energy = ' + str(initialTotalEnergy))
        print('current time step = ' + str(timeStep))
        currentTotalEnergy = get_total_energy_of_pendulum_mp(origin,
                                                             point1Angle, point2Angle,
                                                             point1AngularVelocity, point2AngularVelocity,
                                                             pendulum1Length, pendulum2Length,
                                                             point1Mass, point2Mass,
                                                             gravity)

        print('current total energy = ' + str(currentTotalEnergy))
        print('energy error = ' + str(abs(currentTotalEnergy - initialTotalEnergy)))
        print('')
        timeSinceLastReport = time.time()

# Print the results.
finalTotalEnergy = get_total_energy_of_pendulum_mp(origin,
                                                   point1Angle, point2Angle,
                                                   point1AngularVelocity, point2AngularVelocity,
                                                   pendulum1Length, pendulum2Length,
                                                   point1Mass, point2Mass,
                                                   gravity)

print('final energy error = ' + str(abs(initialTotalEnergy - finalTotalEnergy)))
print('elapsed simulation time = ' + str(elapsedTime))
print('computation time seconds = ' + str(time.time() - start))

max_stars_repo_name: tryabin/double-pendulum-fractal

max_issues_repo_name: tryabin/double-pendulum-fractal

max_forks_repo_name: tryabin/double-pendulum-fractal

reason: from mpmath


text:
__all__ = ['train', 'test', 'classify', 'plot_latent', 'plot_test']

import numpy as np
import matplotlib.pyplot as plt
from copy import deepcopy

import torch
import torch as th
import torch.optim as optim
import torch.nn.functional as F
from torch import nn
from torch.utils.tensorboard import SummaryWriter
from torch.utils.data import random_split
from torchvision.datasets import FashionMNIST

from newremagine import vae
from newremagine import replay

def train(
    fraction,
    train_dataset,
    num_episodes=100,
    batch_size=10,
    num_burn=10,
    lr=0.001,
    device="cpu",
    perfect=True,
    replay_name="Replay",
    replay_kwargs=None,
    vae_name="VAE",
    vae_kwargs=None,
):
    """Given a fraction and dataset, train to self-supervised model

    Params
    -----
    fraction : 3-tuple
        A set of three probability values, setting the probability of
        sampling new data, replaying old data, or imagining data and
        trainong the data.
    train_dataset : a torch dataset object
        The data to train on
    num_epsidoes : int
        The fixed number of training trials
    batch_size : int
        The size of batches to use when traning the network
    num_burn : int (> 0)
        The number of episodes before the we try and replar or
        imagine data. Both these need a min number of experiences
        before they could be useful.
    lr : float (> 0)
        The learning rate of the network
    device : str
        The device to use for training. Either 'cpu` or 'cuda:0'.
        See torch docs for more on this.
    perfect : bool
        If True, replay uses exact copies of the data. If False
        replay uses reconstructed data from the network itself
    replay_name : str
        Any name of a classe found in the `newremainge.replay` module
        can be used here.
    replay_kwargs : dict
        Keword arguments to be transparently passed the the Replay
        memory
    replay_name : str
        Any name of a classe found in the `newremainge.vae` module
        can be used here.
    replay_kwargs : dict
        Keword arguments to be transparently passed the the VAE
        network
    """

    # -- Init memories
    if replay_kwargs is None:
        replay_kwargs = {}
    if vae_kwargs is None:
        vae_kwargs = {}

    Replay = getattr(replay, replay_name)
    memory = Replay(**replay_kwargs)

    VAE = getattr(vae, vae_name)
    model = VAE(**vae_kwargs).to(device)
    optimizer = optim.SGD(model.parameters(), lr=lr)

    # -- Everything is sane?
    num_episodes = int(num_episodes)
    num_burn = int(num_burn)
    lr = float(lr)
    if not np.isclose(sum(fraction), 1):
        raise ValueError("fractions must sum to 1")

    # -- !
    options = ["new", "replay", "imagine"]
    batch_idx = 0
    for n in range(num_episodes):
        # Burn in the model, or consider other options?
        if n < num_burn:
            option = "new"
        else:
            option = np.random.choice(options, p=fraction)

        # Make/get training data
        if option == "new":
            # Get new data
            train_batch = [
                train_dataset[i][0]
                for i in range(batch_idx, batch_idx + batch_size)
            ]
            train_batch = torch.stack(train_batch)
            batch_idx += batch_size
            # If the memory is perfect we replay the
            # training data later, but If the memory
            # is imperfect we add the reconstructed data
            # the the memory intead.
            if perfect:
                memory.encode(train_batch)
            else:
                with torch.no_grad():
                    recon_batch, _, _ = model(train_batch)
                    memory.encode(recon_batch)
        elif option == "replay":
            train_batch = memory.sample(1)[0]
        elif option == "imagine":
            # Sampling from the VAE model is a kind of
            # imagination, or so we imagine in here
            train_batch = model.sample(batch_size, device=device)
        else:
            raise ValueError("Invalid learning option")

        # Train the vae/model
        loss = vae.train(train_batch, model, optimizer, device,
                         model.input_dim)

    return model, memory, float(loss)

def test(model, test_dataset, device="cpu"):
    """Test a pre-trained model on a new dataset.

    Params
    -----
    model : torch nn.Module instance
        The model we want to test.
    test_dataset : a torch dataset object
        The data to train on
     device : str
        The device to use for training. Either 'cpu` or 'cuda:0'.
        See torch docs for more on this
    """
    return vae.test(test_dataset, model, device, model.input_dim)

class _Linear(nn.Module):
    def __init__(self, input_dim=20, output_dim=10):
        super(_Linear, self).__init__()
        # Set dims
        self.input_dim = int(input_dim)
        self.output_dim = int(output_dim)
        # Init layers
        self.fc1 = nn.Linear(self.input_dim, self.output_dim)
        self.logprob = nn.LogSoftmax(dim=1)

    def forward(self, x):
        """Classify c"""
        x = self.fc1(x.view(-1, self.input_dim))
        return self.logprob(x)

def classify(model,
             test_dataset,
             latent_dim=20,
             num_episodes=100,
             batch_size=8,
             lr=0.001,
             device="cpu"):
    """Use a latent space to learn a linear classifier

    Param
    -----
    model : torch nn.Module instance
        The model we want to test.
    test_dataset : a torch dataset object
        The data to train on
    latent_dim : int
        The dim of the latent memory in model
    num_epsidoes : int
        The fixed number of training trials
    batch_size : int
        Training batch size
    lr : float >0
        Learning rate
     device : str
        The device to use for training. Either 'cpu` or 'cuda:0'.
        See torch docs for more on this

    Return
    ------
    loss : float
        Total loss on the test partition
    accuracy ; float
        Final accuracy
    """

    # -- Init
    # Pick some random data, for num_episodes
    n = len(test_dataset)
    idx = np.random.randint(0, 5, size=num_episodes * 2)
    class_dataset = torch.utils.data.Subset(test_dataset, indices=idx)

    # Split it in half
    train, test = random_split(class_dataset, [num_episodes, num_episodes])
    train = torch.utils.data.DataLoader(train, batch_size=batch_size)
    test = torch.utils.data.DataLoader(test, batch_size=batch_size)

    # Init the classifier
    linear = _Linear(input_dim=latent_dim * 2,
                     output_dim=len(test_dataset.classes))
    optimizer = optim.SGD(linear.parameters(), lr=lr)

    # -- !
    for data, labels in train:
        # Get latent encode
        with torch.no_grad():
            data = data.to(device)
            z_mu, z_var = model.encode(data)
            data_z = torch.cat([z_mu, z_var], axis=1)

        # Use it to learn a linear model
        probs = linear(data_z)
        loss = F.nll_loss(probs, labels)
        loss.backward()
        optimizer.step()

    # -- Test final accuracy
    test_loss = 0
    test_correct = 0
    total = 0
    for data, labels in test:
        with torch.no_grad():
            # Get latent encode
            data = data.to(device)
            z_mu, z_var = model.encode(data)
            data_z = torch.cat([z_mu, z_var], axis=1)

            # Test it
            probs = linear(data_z)
            test_loss += F.nll_loss(probs, labels)
            total += labels.size(0)
            predicted = torch.argmax(probs, axis=1)
            test_correct += (predicted == labels).sum().item()

    # --
    return linear, test_loss.item(), test_correct / total

def plot_latent(model, n, img_size=28):
    """Display a grid of samples from the latent space.

    Params
    -----
    model : torch nn.Module instance
        The model we want to test.
    n : int
        The size of the grid
    img_size : int
        The size of images in the orginal data (We assume they
        are greyscale)
    """

    # Make a display grid
    figure = np.zeros((img_size * n, img_size * n))

    # MAke data
    x = model.sample(n**2)
    x = x.flatten().numpy().reshape(n**2, img_size, img_size)
    imgs = [x[i, ::] for i in range(n**2)]

    # !
    k = 0
    for i in range(n):
        for j in range(n):
            figure[i * img_size:(i + 1) * img_size,
                   j * img_size:(j + 1) * img_size] = imgs[k]
            k += 1

    plt.figure(figsize=(10, 10))
    plt.imshow(figure, cmap="gray")
    plt.axis("off")
    plt.show()

def plot_test(model, test_dataset, n, img_size=28):
    """Display a grid random samples from the test data.

    Params
    -----
    model : torch nn.Module instance
        The model we want to test.
    test_dataset : a torch dataset object
        The data to train on
    n : int
        The size of the grid
    img_size : int
        The size of images in the orginal data (We assume they
        are greyscale)
    """

    # Make a display grid
    figure = np.zeros((img_size * n, img_size * n))

    # Choose data
    idx = np.random.randint(0, len(test_dataset), size=n**2)

    # !
    imgs = []
    for i in idx:
        with torch.no_grad():
            x, _ = test_dataset[i]
            x_reconst, _, _ = model(x)

        img = x_reconst.unsqueeze(0).numpy().flatten().reshape(
            img_size, img_size)
        imgs.append(deepcopy(img))

    k = 0
    for i in range(n):
        for j in range(n):
            figure[i * img_size:(i + 1) * img_size,
                   j * img_size:(j + 1) * img_size] = imgs[k]
            k += 1

    plt.figure(figsize=(10, 10))
    plt.imshow(figure, cmap="gray")
    plt.axis("off")
    plt.show()

max_stars_repo_name: CoAxLab/newremagine

max_issues_repo_name: CoAxLab/newremagine

max_forks_repo_name: CoAxLab/newremagine


text:
#!/usr/bin/env python3
import sys, os
import subprocess
import numpy as np
import pandas as pd
import pickle

mapping = {'A':0, # using integers to represent nucleotide types
        'C':1,
        'T':2,
        'G':3,
        'H':4,
        '-':5,
        '':6,
        }

# get samples that have valid oil measurement
traits = pd.read_csv('data/Soybean_Traits.csv', header=0, low_memory=False)
traits = traits[~traits['Oil'].isna()]
pheo_names = [x.replace(' ','') for x in list(traits['Plant ID'])] # get sampleID from pheotype record

# read SNP sequences from genotype file of above samples
snp_dic = {} # keys are sampleID and values are their SNP sequences

# the following line is a naive solution to read large files as it takes a lot memory and time
# you may want to change it to a faster loading method
print('reading SNP file, this may take a long time...')
SNPdf = pd.read_csv('data/soysnp50k_wm82.a2_41317.txt', sep='\t', header=0)
snp_names = list(SNPdf)[5:] # get sampleID from snp record

sampleIDs = list(set(snp_names) & set(pheo_names)) # samples that have both pheotype record and oil value

# dictionary: keys are sampleID and values are cooresponding oil production
oil_dic = {}
for index, row in traits.iterrows():
    sampleID=row['Plant ID'].replace(' ','') # formatting ID to match the style in SNP file
    traits.loc[index, 'Plant ID'] = sampleID # change that in dataframe, too
    if sampleID not in sampleIDs:
        continue
    oil = str(row['Oil']) # read oil measurement
    if ';' in oil:  # use mean if multiple readings are provided
        tmp = [float(x) for x in oil.split(';')]
        oil = np.mean(tmp)
    oil_dic[sampleID] = float(oil) # save to dictionary

traits = traits[traits['Plant ID'].isin(sampleIDs)] # filter out samples from dataframe not presented in SNP records

for sampleID in sampleIDs:
    snp_dic[sampleID] = np.asarray([mapping[x] for x in SNPdf[sampleID]], dtype=np.int8) # using int8 saves more space

assert len(sampleIDs) == traits.shape[0] == len(oil_dic) == len(snp_dic)

preserved_col_names = list(SNPdf)[:5] + sampleIDs
SNPdf = SNPdf[preserved_col_names]
SNPdf.to_csv('data/soysnp50k_wm82.a2_41317_Cleaned.txt', sep='\t', index=False)
print('cleaned snp records saved')

traits.to_csv('data/Soybean_Traits_Cleaned.csv', index=False)
print('cleaned phenotype saved')

pickle.dump(snp_dic, open('data/snp.pickle','wb')) # save snp info to a dictionary
print('snp pickle saved')

pickle.dump(oil_dic, open('data/oil.pickle','wb')) # save oil info
print('oil pickle saved')

max_stars_repo_name: mizaimao/soy

max_issues_repo_name: mizaimao/soy

max_forks_repo_name: mizaimao/soy


text:
# network architecture for face recognition is based on ResNet-34.
# the network is trained by Davis King on LFW having 99.38% accuracy.
import face_recognition
import os
import glob
import argparse
import cv2
import csv
import pandas as pd
import numpy as np

# command line argument
parser = argparse.ArgumentParser()
parser.add_argument('-m', '--method', type=str, choices = ['hog', 'cnn'], default = 'hog',
                    help = 'method to use to detect faces in image')

args = parser.parse_args()

print("[+] Using method : " + args.method)

# path of person images
images_path = os.getcwd() + '/images/'
persons = os.listdir(images_path)

# list to save name and encodings
encodings = []
names = []

for person in persons :
    for person_img in glob.glob(images_path + person + '/' + '*.jpg') :
        img_name = os.path.basename(person_img).split('.')[0]

        print("[+] Creating face encoding for " + person + ' and  image : ' + img_name)
        img_read = cv2.imread(person_img)
        img_rgb = cv2.cvtColor(img_read, cv2.COLOR_BGR2RGB)

        # find the dimensions of box cintaining faces
        faces = face_recognition.face_locations(img_rgb,2,model = args.method)

        # check if face found
        if len(faces) != 0 :
            # create face encodings
            encoding = face_recognition.face_encodings(img_rgb,faces,2)
            encodings.append(np.array(*encoding).tolist())
            names.append(person)

        # if faces are not found
        else :
            print("[-] Could not find face for " + person + 'and image : ' + img_name)

# if method is HOG
if args.method == 'hog' :
    with open(os.getcwd() + '/encodings-hog.csv',mode = 'a') as encode_db :
        writer = csv.writer(encode_db)
        writer.writerow(['Name', 'Encodings'])

        for n,e in zip(names,encodings) :
            writer.writerow([n,e])

# if method is cnn
if args.method == 'cnn' :
    with open(os.getcwd() + '/encodings-cnn.csv',mode = 'a') as encode_db :
        writer = csv.writer(encode_db)
        writer.writerow(['Name', 'Encodings'])

        for n,e in zip(names,encodings) :
            writer.writerow([n,e])


text:
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Mon Apr  6 14:09:05 2020

@author: yannis
"""

import torch
import random
from pdb import set_trace as bp
from a2c_ppo_acktr.envs import make_vec_envs
from a2c_ppo_acktr.utils import  get_vec_normalize
import motion_imitation
import time
import numpy as np

def testPolicy(path,scales=None,pol_scales=None):

    processes = 1

    render = True

    seed = 1

    torch.manual_seed(seed)
    np.random.seed(seed)
    random.seed(seed)

    env = make_vec_envs(
        'A1GymEnv-v1',
        seed,
        processes,
        None,
        None,
        device='cpu',
        allow_early_resets=True, render=render)

    env_core = env.venv.venv.envs[0].env.env
    actor_critic,  ob_rms = torch.load(path,map_location=torch.device('cpu'))
    vec_norm = get_vec_normalize(env)
    if vec_norm is not None:
        vec_norm.eval()
        vec_norm.ob_rms = ob_rms
    recurrent_hidden_states = torch.zeros(1,actor_critic.recurrent_hidden_state_size)
    masks = torch.zeros(1, processes)
    #env_core = env.venv.venv.envs[0]
    if processes==1:
        N_sim = 100
        Reward = np.zeros((N_sim,))
        input('press enter')
        n=0
        R=0
        obs=env.reset()
        while n<N_sim:
            if pol_scales is not None:
                obs[:,-4:] = torch.FloatTensor(pol_scales)
            with torch.no_grad():
                value, action, _, recurrent_hidden_states = actor_critic.act(obs,recurrent_hidden_states,masks, deterministic = True )
            obs, reward, done, _ = env.step(action[0])
            if pol_scales is not None:
                obs[:,-4:] = torch.FloatTensor(pol_scales)
            #env_core.cam_track_torso_link()
            R+=reward
            #control_steps +=1
            time.sleep(5*1.0/240.0)
            if done:
                n+=1
                Reward[n]=R
                print('Reward: ',R)
                R=0
                #obs=env.reset()
                #obs[:,-4:] = torch.FloatTensor(pol_scales)
                #input('press enter')

            masks.fill_(0.0 if done else 1.0)
        #print('Scale: ', Scale[j,:], ', total reward:' , Reward)
        input('press enter')
    else:
        N_sim = processes
        TotalReward = np.zeros((processes,))
        obs=env.reset()
        #bp()
        n = 0
        while n<N_sim:
            if pol_scales is not None:
                obs[:,-4:] = torch.FloatTensor(pol_scales) # replace scale in the input of the policy
            with torch.no_grad():
                    value, action, _, recurrent_hidden_states = actor_critic.act(
                    obs, recurrent_hidden_states, masks, deterministic=True)

            obs, reward, done, _ = env.step(action)
            if pol_scales is not None:
                obs[:,-4:] = torch.FloatTensor(pol_scales) # replace scale in the input of the policy

            TotalReward += reward.numpy().flatten()
            for D in done:
                if D:
                    #print(done)
                    n+=1
            masks = torch.FloatTensor(
                    [[0.0] if done_ else [1.0] for done_ in done])
        print('TotalReward: ', TotalReward, flush=True)
        AverageTotalReward = np.mean(TotalReward)
        Std = np.std(TotalReward)
        #print(TotalReward)
        print('Av. Total reward: ',AverageTotalReward, ', std: ',Std,', virtual scale: ', obs[0,-4:], flush=True)

    #bp()

        N_sim = processes
        TotalReward = np.zeros((processes,))
        obs=env.reset()
        #bp()
        n = 0
        while n<N_sim:
            if pol_scales is not None:
                obs[:,-4:] = torch.FloatTensor(pol_scales) # replace scale in the input of the policy
            with torch.no_grad():
                    value, action, _, recurrent_hidden_states = actor_critic.act(
                    obs, recurrent_hidden_states, masks, deterministic=True)

            obs, reward, done, _ = env.step(action)
            if pol_scales is not None:
                obs[:,-4:] = torch.FloatTensor(pol_scales) # replace scale in the input of the policy

            TotalReward += reward.numpy().flatten()
            for D in done:
                if D:
                    #print(done)
                    n+=1
            masks = torch.FloatTensor(
                    [[0.0] if done_ else [1.0] for done_ in done])
        print('TotalReward: ', TotalReward, flush=True)
        AverageTotalReward = np.mean(TotalReward)
        Std = np.std(TotalReward)
        #print(TotalReward)
        print('Av. Total reward: ',AverageTotalReward, ', std: ',Std,', virtual scale: ', obs[0,-4:], flush=True)
    env.close()
#bp()

if __name__ == '__main__':
    scales = None
    pol_scales = None
    #path = '/home/yannis/Repositories/motion_imitation/12_03_nominal_policy/ppo/A1GymEnv-v1.pt'
    #path = '/home/yannis/Repositories/motion_imitation/12_11_nominal_policy/ppo/A1GymEnv-v1.pt'
    path = '/home/yannis/Repositories/motion_imitation/12_18_nominal_policy/ppo/A1GymEnv-v1.pt'
    testPolicy(path,scales,pol_scales)


text:
import os
import cv2
from skimage import io
import tensorflow as tf
import sys
sys.path.append('../')
import numpy as np
import matplotlib.pyplot as plt
import argparse
import segmentation_models_v1 as sm
from segmentation_models_v1 import Unet, Linknet, PSPNet, FPN, AtUnet
from unet import unet as unet_std

from helper_function import plot_history_flu, plot_set_prediction, generate_folder, plot_prediction_zx, plot_prediction_live
from helper_function import precision, recall, f1_score, calculate_psnr, calculate_pearsonr
from helper_function import plot_flu_prediction, plot_psnr_histogram

sm.set_framework('tf.keras')
import glob
from natsort import natsorted
import time

def read_txt(txt_dir):
    lines = []
    with open(txt_dir, 'r+') as f:
        lines = [line.strip() for line in f.readlines()]
    return lines

parser = argparse.ArgumentParser()
parser.add_argument("--model_file", type=str, default = 'live_models.txt')
parser.add_argument("--index", type=int, default = 0)
parser.add_argument("--gpu", type=str, default = '2')
parser.add_argument("--ch_in", type=int, default = 3)
args = parser.parse_args()
print(args)

os.environ["CUDA_VISIBLE_DEVICES"] = args.gpu
model_names, epoch_list = [], []
with open(args.model_file, 'r+') as f:
		lines = f.readlines()
		for i, line in enumerate(lines):
				model_names.append(line.strip().split(',')[0])
				epoch_list.append(int(line.strip().split(',')[1]))

index = args.index
model_name = model_names[index]
epoch = epoch_list[index]

print(model_name)

## parse model name
splits = model_name.split('-')
dataset = 'neuron_wbx1'
val_dim = 608

scale = 100
subset = 'train2'
chi, cho = 3,3

fl_ch = 'fl2'
for v in range(len(splits)):
	if splits[v] == 'net':
		net_arch = splits[v+1]
	elif splits[v] == 'bone':
		backbone = splits[v+1]
	elif splits[v] == 'scale':
		scale = float(splits[v+1])
	elif splits[v] == 'act':
		act_fun = splits[v+1]
	elif splits[v] == 'scale':
		scale = int(splits[v+1])
	elif splits[v] == 'set':
		dataset = splits[v+1]
	elif splits[v] == 'subset':
		subset = splits[v+1]
	elif splits[v] == 'chi':
		chi = int(splits[v+1])
	elif splits[v] == 'cho':
		cho = int(splits[v+1])
	elif splits[v] == 'chf':
		fl_ch = splits[v+1]

model_folder = '/data/2d_models/{}/{}/'.format(dataset, model_name)
chi = args.ch_in if args.ch_in == 1 else chi

dataset = 'live-neuron_float'
DATA_DIR = '/data/datasets/{}/data/'.format(dataset)
if dataset == 'live-neuron_float':
    val_dim = 1760  # 1744
    offset = 8

volume_fns = [fn for fn in os.listdir(DATA_DIR) if 'output' in fn]

## volum formulation
def extract_vol(vol):
	vol_extr = []
	for i in range(vol.shape[0]):
		if i == 0:
			vol_extr.append(vol[i,:,:,0])
			vol_extr.append(vol[i,:,:,1])
		elif i>0 and i< vol.shape[0]-1:
			vol_extr.append(vol[i,:,:,1])
		else:
			vol_extr.append(vol[i,:,:,1])
			vol_extr.append(vol[i,:,:,2])
	return np.stack(vol_extr)

## metrics
SMOOTH= 1e-6; seed = 0
def calculate_pearsonr(vol1, vol2):
	from scipy import stats
	import numpy as np
	SMOOTH_ = np.random.random(vol1.flatten().shape)*SMOOTH
	flat1 = vol1.flatten(); flat2 = vol2.flatten().flatten()
	flat1 = SMOOTH_ + flat1; flat2 = SMOOTH_ + flat2
	score, _= stats.pearsonr(flat1, flat2)
	return score

def calculate_psnr(vol1, vol2):
	print('value in map: max {}, min {}'.format(vol1.max(),vol2.min()))
	mse = np.mean((vol1-vol2)**2)
	psnr_score = 20*np.log10(255/np.sqrt(mse))
	return psnr_score

# classes for data loading and preprocessing
class Dataset:
    """CamVid Dataset. Read images, apply augmentation and preprocessing transformations.

    Args:
        images_dir (str): path to images folder
        masks_dir (str): path to segmentation masks folder
        class_values (list): values of classes to extract from segmentation mask
        augmentation (albumentations.Compose): data transfromation pipeline
            (e.g. flip, scale, etc.)
        preprocessing (albumentations.Compose): data preprocessing
            (e.g. noralization, shape manipulation, etc.)

    """

    def __init__(
            self,
            images_dir,
            masks_dir1,
            masks_dir2,
            fl_ch = None,
            scale = 1.0,
            channels = [3,3],
            nb_data=None,
            augmentation=None,
            preprocessing=None,
    ):
        #self.ids = os.listdir(images_dir)
        id_list = natsorted(os.listdir(images_dir))
        if nb_data ==None:
            self.ids = id_list
        else:
            self.ids = id_list[:int(min(nb_data,len(id_list)))]
        self.images_fps = [os.path.join(images_dir, image_id) for image_id in self.ids]
        self.masks1_fps = [os.path.join(masks_dir1, image_id) for image_id in self.ids]
        self.masks2_fps = [os.path.join(masks_dir2, image_id) for image_id in self.ids]
        print('Load files: image {}, fl1: {}, fl2:{}'.format(len(self.images_fps),len(self.masks1_fps),len(self.masks2_fps)))
        self.scale = scale
        self.augmentation = augmentation
        self.preprocessing = preprocessing
        self.channels = channels
        self.fl_ch = fl_ch

    def __getitem__(self, i):

        # load image and fl1 or fl2 or both
#         image = io.imread(self.images_fps[i])
#         if self.fl_ch == 'fl1':
#             mask  = io.imread(self.masks1_fps[i])
#             mask = mask/255.*self.scale
#         elif self.fl_ch == 'fl2':
#             mask  = io.imread(self.masks2_fps[i])
#             mask = mask/255.*self.scale
#         elif self.fl_ch == 'fl12':
#             mask1  = io.imread(self.masks1_fps[i])
#             mask2  = io.imread(self.masks2_fps[i])
#             mask = np.stack([mask1[:,:,1], mask2[:,:,1]], axis = -1)
#             mask = mask/255.*self.scale
        # load image and fl1 or fl2 or both
        image = np.load(self.images_fps[i]) * 255.
        if self.fl_ch == 'fl1':
            mask  = np.load(self.masks1_fps[i])
            mask = mask *self.scale
        elif self.fl_ch == 'fl2':
            mask  = np.load(self.masks2_fps[i])
            mask = mask *self.scale
        elif self.fl_ch == 'fl12':
            mask1  = np.load(self.masks1_fps[i])
            mask2  = np.load(self.masks2_fps[i])
            mask = np.stack([mask1[:,:,1], mask2[:,:,1]], axis = -1)
            mask = mask *self.scale

        # decide the input and output channels
        if self.channels[0] == 1:
            image[:,:,0], image[:,:,2] = image[:,:,1], image[:,:,1]
        elif self.channels[0] == 2:
            image[:,:,2] = image[:,:,1]

        if self.channels[1] == 1 and not (self.fl_ch=='fl12'):
            mask = mask[:,:,1:2]

        # apply augmentations
        if self.augmentation:
            sample = self.augmentation(image=image, mask=mask)
            image, mask = sample['image'], sample['mask']

        # apply preprocessing
        if self.preprocessing:
            sample = self.preprocessing(image=image, mask=mask)
            image, mask = sample['image'], sample['mask']

        return image, mask

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

class Dataloder(tf.keras.utils.Sequence):
    """Load data from dataset and form batches

    Args:
        dataset: instance of Dataset class for image loading and preprocessing.
        batch_size: Integet number of images in batch.
        shuffle: Boolean, if `True` shuffle image indexes each epoch.
    """

    def __init__(self, dataset, batch_size=1, shuffle=False):
        self.dataset = dataset
        self.batch_size = batch_size
        self.shuffle = shuffle
        self.indexes = np.arange(len(dataset))

        self.on_epoch_end()

    def __getitem__(self, i):

        # collect batch data
        start = i * self.batch_size
        stop = (i + 1) * self.batch_size
        data = []
        for j in range(start, stop):
            data.append(self.dataset[j])

        # transpose list of lists
        batch = [np.stack(samples, axis=0) for samples in zip(*data)]

        return (batch[0], batch[1])

    def __len__(self):
        """Denotes the number of batches per epoch"""
        return len(self.indexes) // self.batch_size

    def on_epoch_end(self):
        """Callback function to shuffle indexes each epoch"""
        if self.shuffle:
            self.indexes = np.random.permutation(self.indexes)

import albumentations as A

# define heavy augmentations
def get_training_augmentation(dim, rot = 0):
    train_transform = [
        A.PadIfNeeded(min_height=dim, min_width=dim, always_apply=True, border_mode=0),
        A.RandomCrop(height=dim, width=dim, always_apply=True),
    ]
    return A.Compose(train_transform)

def get_validation_augmentation(dim = 992):
    """Add paddings to make image shape divisible by 32"""
    test_transform = [
        A.PadIfNeeded(dim, dim)
    ]
    return A.Compose(test_transform)

def get_preprocessing(preprocessing_fn):
    """Construct preprocessing transform
    Args:
        preprocessing_fn (callbale): data normalization function
            (can be specific for each pretrained neural network)
    Return:
        transform: albumentations.Compose
    """

    _transform = [
        A.Lambda(image=preprocessing_fn),
    ]
    return A.Compose(_transform)

# network
best_weight = model_folder+'/ready_model.h5'
preprocess_input = sm.get_preprocessing(backbone)

#create model
n_classes = cho if not fl_ch == 'fl12' else 2
activation = 'relu'
net_func = globals()[net_arch]
model = net_func(backbone, classes=n_classes, encoder_weights=None, activation=activation)

#load best weights
model.load_weights(best_weight)

## save model
# model.save(model_folder+'/ready_model.h5')

subset = 'live-neuron'
vol_fns = volume_fns
for vol_fn in vol_fns:
		if not os.path.exists(os.path.join(DATA_DIR, vol_fn)):
				continue
		print('{}: {}'.format(subset, vol_fn))
		X_dir = os.path.join(DATA_DIR, vol_fn, 'phase')
		Y1_dir = os.path.join(DATA_DIR, vol_fn, 'phase')
		Y2_dir = os.path.join(DATA_DIR, vol_fn, 'phase')

		test_dataset = Dataset(
			X_dir,
			Y1_dir,
			Y2_dir,
			fl_ch = fl_ch,
			channels = [chi, cho],
			scale = scale,
			augmentation=get_validation_augmentation(val_dim),
			preprocessing=get_preprocessing(preprocess_input),
		)

		print(test_dataset[0][0].shape, test_dataset[0][1].shape)

		test_dataloader = Dataloder(test_dataset, batch_size=1, shuffle=False)

		# prediction and ground truth
		pr_masks = model.predict(test_dataloader)
		print(pr_masks.shape)
		images = []
		for i in range(len(test_dataset)):
			image, _ = test_dataset[i];images.append(image)
		images = np.stack(images)

		# scale to [0,255]
# 		gt_masks = np.uint8(gt_masks/scale*255); gt_masks = gt_masks.squeeze()
		pr_masks = pr_masks/scale*255.; pr_masks = pr_masks.squeeze()
		images = np.uint8(images)

		# obtain the images, GT, and prediction
		ph_vol = extract_vol(images)
		ph_vol = ph_vol[:,offset:-offset,offset:-offset]

		# psnr and pearsonr correlation
		if fl_ch == 'fl1':
				if cho == 3:
						pr_vol = extract_vol(pr_masks)
				else:
						pr_vol = pr_masks
		elif fl_ch == 'fl2':
				if cho == 3:
						pr_vol2 = extract_vol(pr_masks)
				else:
						pr_vol2 = pr_masks
		elif fl_ch == 'fl12':
				pr_vol = pr_masks[:,:,:,0]
				pr_vol2 = pr_masks[:,:,:,1]

		if fl_ch == 'fl12' or fl_ch == 'fl1':
				pr_vol = pr_vol[:,offset:-offset,offset:-offset]
		if fl_ch == 'fl12' or fl_ch == 'fl2':
				pr_vol2 = pr_vol2[:,offset:-offset,offset:-offset]

		# save prediction
		pred_save = True
		if pred_save:
				pr_vol_dir = model_folder+'/pred_fl1_fl2_live' if chi == 3 else model_folder+'/pred_fl1_fl2_live_one_input'
				generate_folder(pr_vol_dir)
				if fl_ch == 'fl12' or fl_ch == 'fl1':
						np.save(os.path.join(pr_vol_dir,'Pr1_{}.npy'.format(vol_fn)), pr_vol)
						print('FL1: {}'.format(pr_vol.shape))
				if fl_ch == 'fl12' or fl_ch == 'fl2':
						np.save(os.path.join(pr_vol_dir,'Pr2_{}.npy'.format(vol_fn)), pr_vol2)
						print('FL2: {}'.format(pr_vol2.shape))

		# save prediction examples
		prediction_dir = model_folder+'/pred_examples_live'
		generate_folder(prediction_dir)
		plot_fig_file = prediction_dir+'/{}_fl1.png'.format(vol_fn)
		plot_fig_file2 = prediction_dir+'/{}_fl2.png'.format(vol_fn)
		nb_images = 10
		if fl_ch == 'fl12' or fl_ch == 'fl1':
				z_index = int(pr_vol.shape[0]/2)
# 				plot_flu_prediction(plot_fig_file, ph_vol, pr_vol, pr_vol, nb_images)
				plot_prediction_live(plot_fig_file, ph_vol, pr_vol, pr_vol, z_index)
		if fl_ch == 'fl12' or fl_ch == 'fl2':
				z_index = int(pr_vol2.shape[0]/2)
				plot_prediction_live(plot_fig_file2, ph_vol, pr_vol2, pr_vol2, z_index)
# 				plot_flu_prediction(plot_fig_file2, ph_vol, pr_vol2, pr_vol2, nb_images)

reason: import numpy,from scipy


text:
def drhankelapprox(H,r,solver = None, gamma = 1,rho = 1,Z0 = None,tol = None):
    """
    Douglas-Rachford proximal splitting for low-rank
    approximation with Hankel constraint through the low-rank inducing
    Frobenius norm and non-convex Douglas-Rachford.

    M,rankM,e,D,Z_fix,iter = DRHANKELAPPROX(H,r) determines a Frobenius norm
    low-rank Hankel approximation M through the Frobenius norm low-rank
    inducing norm and Douglas-Rachford splitting, i.e., M is the solution
    to

    minimize 0.5*||H||_{ell_2}^2 - trace(M'H) + 0.5*||M||_{ell_2,r*}^2
            s.t. M is Hankel

    Further,
        1. rankM = rank(M)
        2. err = norm(H-M,'fro')
        3. D = solution of the dual problem
        4. Z_fix = fix point of the Douglas-Rachford iterations
        5. iter = total number of Douglas-Rachford iterations

    ... = DRHANKELAPPROX(H,r,option) allows to specify further options:
        1. ... = DRHANKELAPPROX(H,r,...,solver = 'NDR',...) changes solver
        to use the non-convex Douglas-Rachford.
        2. ... = DRHANKELAPPROX(H,r,...,gamma = gamma_val,...) multiplies the
        objective functions with gamma_val when determining the prox of them.
        Default value: gamma_val = 1.
        3. ... = DRHANKELAPPROX(H,r,...,rho = rho_val,...) set the step length
        update of the fix-point update, i.e.
        Z_{k+1} = Z_k + rho_val*(Y_k - X_k), where 0 < rho_val < 2.
        Default value: rho_val = 1.
        4. ... = DRHANKELAPPROX(H,r,...,Z0 = Z0_val,...) sets the initial value of
        the fix-point iteration, i.e., Z_{0} = Z0_val.
        Default choice: Z0 = zeros(H.shape).
        5. ... = DRHANKELAPPROX(H,r,...,tol = tol_val,...) sets the relative
        tolerance of:
            + The numerical rank: rankM = rank(M/norm(H,'fro'),tol_val)
            + Iterations stop: Stop if (Y_k -X_k)/norm(H,'fro') < tol_val
            + Zero values: E.g., D(abs(D/norm(H,'fro')) < tol_val) = 0
        Default value: tol_val = sqrt(eps), where eps denotes the machine precision constant.

    ---------------
    References:
    - C. Grussler and A. Rantzer and P. Giselsson (2018):
    "Low-Rank Optimization with Convex Constraints",
    IEEE Transactions on Automatic Control, DOI: 10.1109/TAC.2018.2813009.

    - C. Grussler and P. Giselsson (2018):
    "Low-Rank Inducing Norms With Optimality Interpreations",
    SIAM J. Optim., 28(4), pp. 3057–3078.

    - C. Grussler and P. Giselsson (2017):
    "Local convergence of proximal splitting methods for rank constrained
    problems", pp. 702-708, IEEE 56th Annual Conference on Decision and Control
    (CDC), DOI: 10.1109/CDC.2017.8263743.

    - C. Grussler (2017):
    "Rank reduction with convex constraints", PhD Thesis,
    Department of Automatic Control, Lund Institute of Technology,
    Lund University, ISBN 978-91-7753-081-7.
    """

    import numpy as np
    from functools import partial
    from lripy.projhankel import projhankel
    from lripy.proxnonconv_square import proxnonconv_square
    from lripy.proxnormrast_square import proxnormrast_square
    from lripy.dr import dr

    # Define default tolerance value
    if tol is None:
        tol = np.sqrt(np.finfo(float).eps)

    dim = H.shape

    if r > min(dim):
        ValueError('r is larger than min(H.shape)')

    # Define absolute tolerance
    tol_H = np.linalg.norm(H)*tol

    ## Start Douglas-Rachford iterations
    # Choose between convex and non-convex Douglas-Rachford
    if  solver == 'NDR':
        M,Z_fix,iter,D = dr(partial(proxnonconv_square,r=r,p=2),partial(projhankel,H = -H),dim,tol =tol_H,gamma = gamma,Z0 = Z0,rho = rho)
    else:
        M,Z_fix,iter,D = dr(partial(proxnormrast_square,r=r,p=2),partial(projhankel,H = -H),dim,tol =tol_H,gamma = gamma,Z0 = Z0,rho = rho)

    # Compute error and rank of the approximation
    err = np.linalg.norm(M-H)
    rankM = np.linalg.matrix_rank(M,tol_H)

    return M,rankM,err,D,Z_fix,iter

max_stars_repo_name: LowRankOpt/LRINormPy

max_issues_repo_name: LowRankOpt/LRINormPy

max_forks_repo_name: LowRankOpt/LRINormPy


text:
import pandas as pd
import config as cfg
import requests
import pymysql
import numpy as np
import csv

api_key = cfg.settings["census"]["key"]

connection = cfg.get_db_connection()

try:
    with connection.cursor() as cursor:
        # Read a single record
        cursor.execute("SELECT Geoid FROM community_units WHERE community_id=1")
        rows = cursor.fetchall()
        print(rows)
finally:
    connection.close()

list_of_bgs=[]
for i in rows:
    i = ''.join(filter(lambda x: x in '.0123456789', str(i)))
    i = str(i)
    list_of_bgs.append(i)

# list_of_bgs = ['371830501001','371830507001','371830507002']
#list_of_bgs = ['371830501001']

def pull_census(geoid):
    state = geoid[0:2]
    county = geoid[2:5]
    tract = geoid[5:11]
    blockgroup = geoid[11:]

    url = ("https://api.census.gov/data/2010/sf1?get=NAME,P0050001,P0050002,P0050003,P0050004,P0050005,"+
            "P0050006,P0050007,P0050008,P0050009,P0050010,P0050011,P0050012,P0050013,P0050014,P0050015,"+
            "P0050016,P0050017,H0110001,H0110002,H0110003,H0110004&for=block%20group:"+ blockgroup +
            "&in=state:" + state + "%20county:" + county + "%20tract:" + tract + "&key=" + api_key)
    html = requests.get(url).json()
    bg_data = pd.DataFrame(html, columns=html[0])
    return bg_data[1:]

census_data = pull_census(list_of_bgs[0])

if len(list_of_bgs) > 1:
    for bg in list_of_bgs[1:]:
        next_bg = pull_census(bg)
        census_data = pd.concat([census_data,next_bg])

        # Variables of interest:
        # - P0050001-P00500018: Race/Ethnicity
        # - H0100001-H0100004: Population in Owned v. Rented units

#asnumeric
census_data['P0050001']=census_data['P0050001'].astype(str).astype(int)
census_data['P0050002']=census_data['P0050002'].astype(str).astype(int)
census_data['P0050003']=census_data['P0050003'].astype(str).astype(int)
census_data['P0050004']=census_data['P0050004'].astype(str).astype(int)
census_data['P0050005']=census_data['P0050005'].astype(str).astype(int)
census_data['P0050006']=census_data['P0050006'].astype(str).astype(int)
census_data['P0050007']=census_data['P0050007'].astype(str).astype(int)
census_data['P0050008']=census_data['P0050008'].astype(str).astype(int)
census_data['P0050009']=census_data['P0050009'].astype(str).astype(int)
census_data['P0050010']=census_data['P0050010'].astype(str).astype(int)
census_data['P0050011']=census_data['P0050011'].astype(str).astype(int)
census_data['P0050012']=census_data['P0050012'].astype(str).astype(int)
census_data['P0050013']=census_data['P0050013'].astype(str).astype(int)
census_data['P0050014']=census_data['P0050014'].astype(str).astype(int)
census_data['P0050015']=census_data['P0050015'].astype(str).astype(int)
census_data['P0050016']=census_data['P0050016'].astype(str).astype(int)
census_data['P0050017']=census_data['P0050017'].astype(str).astype(int)
census_data['H0110001']=census_data['H0110001'].astype(str).astype(int)
census_data['H0110002']=census_data['H0110002'].astype(str).astype(int)
census_data['H0110003']=census_data['H0110003'].astype(str).astype(int)
census_data['H0110004']=census_data['H0110004'].astype(str).astype(int)
census_data['Hispanic']=census_data['P0050001']-census_data['P0050002']

#calcs
pct_white_non_hispanic=census_data['P0050003'].sum()/census_data['P0050001'].sum()
pct_black_non_hispanic=census_data['P0050004'].sum()/census_data['P0050001'].sum()
pct_hispanic=census_data['Hispanic'].sum()/census_data['P0050001'].sum()
pct_AINI_non_hispanic=census_data['P0050005'].sum()/census_data['P0050001'].sum()
pct_asian_non_hispanic=census_data['P0050006'].sum()/census_data['P0050001'].sum()
pct_NHPI_non_hispanic=census_data['P0050007'].sum()/census_data['P0050001'].sum()
pct_other_non_hispanic=census_data['P0050008'].sum()/census_data['P0050001'].sum()
pct_multi_non_hispanic=census_data['P0050009'].sum()/census_data['P0050001'].sum()
pct_renters=census_data['H0110004'].sum()/census_data['H0110001'].sum()

d = {"race:pct_white_non_hispanic":pct_white_non_hispanic, "race:pct_black_non_hispanic": pct_black_non_hispanic, "race:pct_AINI_non_hispanic": pct_AINI_non_hispanic, "race:pct_asian_non_hispanic": pct_asian_non_hispanic,"race:pct_hispanic":pct_hispanic, "race:pct_other_non_hispanic":pct_other_non_hispanic, "race:pct_multi_non_hispanic":pct_multi_non_hispanic, "housing:pct_renters":pct_renters }

(pd.DataFrame.from_dict(data=d, orient='index')
   .to_csv('dict_file.csv', header=False))

max_stars_repo_name: gerrymandr/community-builder

max_issues_repo_name: gerrymandr/community-builder

max_forks_repo_name: gerrymandr/community-builder


text:
import numpy as np
import trimesh
import pandas as pd
from scipy.stats import entropy
from skimage import color
from sklearn.preprocessing import scale
from nss_functions import *
import os

def get_feature_vector(objpath):
  print("Begin loading mesh.")
  mesh = trimesh.load(objpath)

  print("Begin geometry feature extraction.")
  # extract geometry features
  to_origin, extents = trimesh.bounds.oriented_bounds(mesh)
  print("Mean curvature")
  r = max(extents) * 0.003
  #mean_curvature = np.array(trimesh.curvature.discrete_mean_curvature_measure(mesh, mesh.vertices, r))[0]
  mean_curvature = np.array(trimesh.curvature.discrete_gaussian_curvature_measure(mesh, mesh.vertices, r))
  print("Dihedral angle")
  dihedral_angle = mesh.face_adjacency_angles
  print("Face")
  face_area = mesh.area_faces
  face_angle = np.array([i for item in mesh.face_angles for i in item])

  print("Begin color feature extraction.")
  # extract color features
  rgb_color = mesh.visual.vertex_colors/255
  rgb_color = np.delete(rgb_color,3,axis=1)
  lab_color = color.rgb2lab(rgb_color)
  l = lab_color[:,0]
  a = lab_color[:,1]
  b = lab_color[:,2]

  print("Begin NSS parameters estimation.")
  # computer nss parameters
  nss_params = []
  i = 0
  for tmp in [l,a,b,mean_curvature,dihedral_angle,face_angle,face_area]:
      i = i + 1
      print(i)
      params = get_nss_param(tmp)
      #flatten the feature vector
      nss_params = nss_params + [i for item in params for i in item]
  return nss_params

objpath = "Chameleon_SimpMesh_4.obj"
features = get_feature_vector(objpath)
print(features)


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

image_size = cfg.image_Width

file = "C:\\dataset\\VOC2012\\JPEGImages\\2011_006205.jpg"

def image_read(imname, flipped=False):  # 读取图片
    image = cv2.imread(imname)
    original_h = image.shape[0]  # height = row
    original_w = image.shape[1]  # width = column(=col)
    print(original_w)
    image = cv2.resize(image, (image_size, image_size))  # resize大小
    image = cv2.cvtColor(image, cv2.COLOR_BGR2RGB).astype(np.float32)
    image = (image / 255.0) * 2.0 - 1.0
    print(f"image name = {imname}")
    plt.imshow(image)
    plt.show()
    if flipped:
        image = image[:, ::-1, :]

    return image

if __name__ == "__main__":
    image = image_read(file)


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

import numpy as np

import gzip
import bz2

import argparse
import logging
import sys

def iopen(file, *args, **kwargs):
    _open = open
    if file.endswith('.gz'):
        _open = gzip.open
    elif file.endswith('.bz2'):
        _open = bz2.open
    return _open(file, *args, **kwargs)

def read_triples(path):
    triples = []
    with iopen(path, 'rt') as f:
        for line in f.readlines():
            s, p, o = line.split()
            triples += [(s.strip(), p.strip(), o.strip())]
    return triples

def main(argv):
    argparser = argparse.ArgumentParser('Knowledge Graph Splitter',
                                        formatter_class=argparse.ArgumentDefaultsHelpFormatter)

    argparser.add_argument('triples', type=str, default='/dev/stdin')

    argparser.add_argument('--kept', type=argparse.FileType('w'), default='/dev/stdout')
    argparser.add_argument('--held-out', type=argparse.FileType('w'), default='/dev/stdout')
    argparser.add_argument('--held-out-size', action='store', type=int, default=1000)

    argparser.add_argument('--seed', action='store', type=int, default=0)

    args = argparser.parse_args(argv)

    triples_path = args.triples

    kept_fd = args.kept

    held_out_fd = args.held_out
    held_out_size = args.held_out_size

    assert held_out_size > 0

    seed = args.seed

    triples = read_triples(triples_path)
    nb_triples = len(triples)

    nb_kept_triples = nb_triples - held_out_size
    assert nb_kept_triples > 0

    random_state = np.random.RandomState(seed)
    permutation = random_state.permutation(nb_triples).tolist()

    shuffled_triples = [triples[i] for i in permutation]

    kept_triples = shuffled_triples[:nb_kept_triples]
    held_out_triples = shuffled_triples[nb_kept_triples:]

    kept_fd.writelines(['\t'.join(symbol for symbol in triple) + '\n' for triple in kept_triples])
    held_out_fd.writelines(['\t'.join(symbol for symbol in triple) + '\n' for triple in held_out_triples])

if __name__ == '__main__':
    logging.basicConfig(level=logging.DEBUG)
    main(sys.argv[1:])


text:
#!/usr/bin/env python3

from phe import paillier
import numpy as np

# Generate the public and private key used for Paillier encryption and decryption
public_key, private_key = paillier.generate_paillier_keypair()

# List of numbers we want to encrypt
secret_number_list = [3, 5, 7, 9, 11]

# We encrypt the list using the public key
encrypted_number_list = [public_key.encrypt(x) for x in secret_number_list]

def sum_mult(encrypted_number_list, private_key):
    sum = 0
    for i in range(len(encrypted_number_list)) :
        # encrypted_number_list[i] = encrypted_number_list[i] * 2
        sum += encrypted_number_list[i] * 2
    # encrypted_mean = np.sum(encrypted_number_list)

    return private_key.decrypt(sum)

def mean(encrypted_number_list, private_key):
    sum = 0
    for i in range(len(encrypted_number_list)):
        sum += encrypted_number_list[i]
    mean = sum / len(encrypted_number_list)
    return private_key.decrypt(mean)

print("sum of element multiplied by 2 =", sum_mult(encrypted_number_list, private_key))
print("mean of the list is =", mean(encrypted_number_list, private_key))


text:
import numpy as np
import cv2
import os

def distance(v1, v2):
    # Eucledian
    return np.sqrt(((v1-v2)**2).sum())

def knn(train, test, k=5):
    dist = []

    for i in range(train.shape[0]):
        #Get the vector and label
        ix = train[i, :-1]
        iy = train[i, -1]
        # Computhe the distance from test point
        d = distance(test, ix)
        dist.append([d, iy])
    # Sort based on distance and get top k
    dk = sorted(dist, key=lambda x: x[0])[:k]
    # Retrieve only the labels
    labels = np.array(dk)[:, -1]

    # Get frequencies of each label
    output = np.unique(labels, return_counts=True)
    # Find max frequency and corresponding label
    index = np.argmax(output[1])
    return output[0][index]

cap = cv2.VideoCapture(0)
face_cascade = cv2.CascadeClassifier("haarcascade_frontalface_alt.xml") ## Taking the instructions from opencv

dataset_path = "./face_dataset/" ## Path to get the face data

face_data = []
labels = []
class_id = 0 ## Labels for every given file
names = {} ## Mapping between if and name

# Dataset preparation
for fx in os.listdir(dataset_path): ## Loop on the faces data folder.
    if fx.endswith('.npy'):
        names[class_id] = fx[:-4]
        data_item = np.load(dataset_path + fx)
        face_data.append(data_item)

        target = class_id * np.ones((data_item.shape[0],))
        class_id += 1
        labels.append(target)

face_dataset = np.concatenate(face_data, axis=0)
face_labels = np.concatenate(labels, axis=0).reshape((-1, 1))
print(face_labels.shape)
print(face_dataset.shape)

trainset = np.concatenate((face_dataset, face_labels), axis=1)
print(trainset.shape)

font = cv2.FONT_HERSHEY_SIMPLEX

while True:
    ret, frame = cap.read()
    if ret == False:
        continue

    # Convert frame to grayscale
    gray = cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY)

    # Detect multi faces in the image
    faces = face_cascade.detectMultiScale(gray, 1.3, 5)

    for face in faces:
        x, y, w, h = face

        # Get the face RDI
        offset = 5
        face_section = frame[y-offset:y+h+offset, x-offset:x+w+offset]
        face_section = cv2.resize(face_section, (100, 100))

        out = knn(trainset, face_section.flatten())

        # Draw rectangle in hte original image
        cv2.putText(frame, names[int(out)],(x,y-10), cv2.FONT_HERSHEY_SIMPLEX, 1,(255,0,0),2,cv2.LINE_AA)
        cv2.rectangle(frame, (x,y), (x+w,y+h), (255,255,255), 2)

    cv2.imshow("Faces", frame)

    if cv2.waitKey(1) & 0xFF == ord('q'):
        break

cap.release()
cv2.destroyAllWindows() ## Close the camera.