text:
import sys
import os
import re
import time
import timeit
import numpy as np
import numpy.random as rdn
from qtools.qtpy import QtCore, QtGui
from qtools.qtpy.QtCore import Qt, pyqtSignal
from qtools import show_window
from galry import DEBUG, log_debug, log_info, log_warn
try:
    from qtools.qtpy.QtOpenGL import QGLWidget, QGLFormat
except Exception as e:
    log_warn((("The Qt-OpenGL bindings are not available. "
    "On Ubuntu, please install python-qt4-gl (PyQt4) or "
    "python-pyside.qtopengl (PySide). "
    "Original exception was: %s" % e)))
    # mock QGLWidget
    class QGLWidget(QtGui.QWidget):
        def initializeGL(self):
            pass
        def paintGL(self):
            pass
        def updateGL(self):
            pass
        def resizeGL(self):
            pass
    QGLFormat = None
from galry import get_cursor, FpsCounter, PaintManager, \
    InteractionManager, BindingManager, \
    UserActionGenerator, PlotBindings, Bindings, FpsCounter, \
    show_window, get_icon

__all__ = [
'GalryWidget',
'GalryTimerWidget',
'AutodestructibleWindow',
'create_custom_widget',
'create_basic_window',
'show_basic_window',
]

# # DEBUG: raise errors if Numpy arrays are unnecessarily copied
# from OpenGL.arrays import numpymodule
# try:
    # numpymodule.NumpyHandler.ERROR_ON_COPY = True
# except Exception as e:
    # print "WARNING: unable to set the Numpy-OpenGL copy warning"

# Set to True or to a number of milliseconds to have all windows automatically
# killed after a fixed time. It is useful for automatic debugging or
# benchmarking.
AUTODESTRUCT = False
DEFAULT_AUTODESTRUCT = 1000

# Display the FPS or not.
DISPLAY_FPS = DEBUG == True

# Default manager classes.
DEFAULT_MANAGERS = dict(
    paint_manager=PaintManager,
    binding_manager=BindingManager,
    interaction_manager=InteractionManager,
)

# Main Galry class.
class GalryWidget(QGLWidget):
    """Efficient interactive 2D visualization widget.

    This QT widget is based on OpenGL and depends on both PyQT (or PySide)
    and PyOpenGL. It implements low-level mechanisms for interaction processing
    and acts as a glue between the different managers (PaintManager,
    BindingManager, InteractionManager).

    """

    w = 600.
    h = 600.

    # Initialization methods
    # ----------------------
    def __init__(self, format=None, autosave=None, getfocus=True, **kwargs):
        """Constructor. Call `initialize` and initialize the companion classes
        as well."""
        if format is not None:
            super(GalryWidget, self).__init__(format)
        else:
            super(GalryWidget, self).__init__()

        self.initialized = False
        self.just_initialized = False

        self.i = 0

        # background color as a 4-tuple (R,G,B,A)
        self.bgcolor = (0, 0, 0, 0)
        self.autosave = None

        # default window size
        # self.width, self.height = 600, 600

        # FPS counter, used for debugging
        self.fps_counter = FpsCounter()
        self.display_fps = DISPLAY_FPS
        self.activate3D = None

        # widget creation parameters
        self.bindings = None
        self.companion_classes_initialized = False

        # constrain width/height ratio when resizing of zooming
        self.constrain_ratio = False
        self.constrain_navigation = False
        self.momentum = False
        self.activate_help = True
        self.activate_grid = False
        self.block_refresh = False

        # Capture keyboard events.
        if getfocus:
            self.setFocusPolicy(Qt.WheelFocus)

        # Capture mouse events.
        self.setMouseTracking(True)

        # Capture touch events.
        self.setAcceptTouchEvents = True
        self.grabGesture(QtCore.Qt.PinchGesture)
        self.mouse_blocked = False  # True during a pinch gesture

        # Initialize the objects providing the core features of the widget.
        self.user_action_generator = UserActionGenerator()

        self.is_fullscreen = False

        self.events_to_signals = {}

        # keyword arguments without "_manager" => passed to initialize
        self.initialize(**kwargs)

        # initialize companion classes if it has not been done in initialize
        if not self.companion_classes_initialized:
            self.initialize_companion_classes()
        self.initialize_bindings()

        # update rendering options
        self.paint_manager.set_rendering_options(
                        activate3D=self.activate3D,
                        constrain_ratio=self.constrain_ratio,
                        )

        self.autosave = autosave

    def set_bindings(self, *bindings):
        """Set the interaction mode by specifying the binding object.

        Several binding objects can be given for the binding manager, such that
        the first one is the currently active one.

        Arguments:
          * bindings: a list of classes instances deriving from
            Bindings.

        """
        bindings = list(bindings)
        if not bindings:
            bindings = [PlotBindings()]
        # if type(bindings) is not list and type(bindings) is not tuple:
            # bindings = [bindings]
        # if binding is a class, try instanciating it
        for i in xrange(len(bindings)):
            if not isinstance(bindings[i], Bindings):
                bindings[i] = bindings[i]()
        self.bindings = bindings

    def set_companion_classes(self, **kwargs):
        """Set specified companion classes, unspecified ones are set to
        default classes.

        Arguments:
          * **kwargs: the naming convention is: `paint_manager=PaintManager`.
            The key `paint_manager` is the name the manager is accessed from
            this widget and from all other companion classes. The value
            is the name of the class, it should end with `Manager`.

        """
        if not hasattr(self, "companion_classes"):
            self.companion_classes = {}

        self.companion_classes.update(kwargs)

        # default companion classes
        self.companion_classes.update([(k,v) for k,v in \
            DEFAULT_MANAGERS.iteritems() if k not in self.companion_classes])

    def initialize_bindings(self):
        """Initialize the interaction bindings."""
        if self.bindings is None:
            self.set_bindings()
        self.binding_manager.add(*self.bindings)

    def initialize_companion_classes(self):
        """Initialize companion classes."""
        # default companion classes
        if not getattr(self, "companion_classes", None):
            self.set_companion_classes()

        # create the managers
        for key, val in self.companion_classes.iteritems():
            log_debug("Initializing '%s'" % key)
            obj = val(self)
            setattr(self, key, obj)

        # link all managers
        for key, val in self.companion_classes.iteritems():
            for child_key, child_val in self.companion_classes.iteritems():
                # no self-reference
                if child_key == key:
                    continue
                obj = getattr(self, key)
                setattr(obj, child_key, getattr(self, child_key))

        self.interaction_manager.constrain_navigation = self.constrain_navigation
        self.companion_classes_initialized = True

    def initialize(self, **kwargs):
        """Initialize the widget.

        Parameters such as bindings, companion_classes can be
        set here, by overriding this method. If initializations must be done
        after companion classes instanciation, then
        self.initialize_companion_classes can be called here.
        Otherwise, it will be called automatically after initialize().

        """
        pass

    def clear(self):
        """Clear the view."""
        self.paint_manager.reset()

    def reinit(self):
        """Reinitialize OpenGL.

        The clear method should be called before.

        """
        self.initializeGL()
        self.resizeGL(self.w, self.h)
        self.updateGL()

    # OpenGL widget methods
    # ---------------------
    def initializeGL(self):
        """Initialize OpenGL parameters."""
        self.paint_manager.initializeGL()
        self.initialized = True
        self.just_initialized = True

    def paintGL(self):
        """Paint the scene.

        Called as soon as the window needs to be painted (e.g. call to
        `updateGL()`).

        This method calls the `paint_all` method of the PaintManager.

        """
        if self.just_initialized:
            self.process_interaction('Initialize', do_update=False)
        # paint fps
        if self.display_fps:
            self.paint_fps()
        # paint everything
        self.paint_manager.paintGL()
        # compute FPS
        self.fps_counter.tick()
        if self.autosave:
            if '%' in self.autosave:
                autosave = self.autosave % self.i
            else:
                autosave = self.autosave
            self.save_image(autosave, update=False)
        self.just_initialized = False
        self.i += 1

    def paint_fps(self):
        """Display the FPS on the top-left of the screen."""
        self.paint_manager.update_fps(int(self.fps_counter.get_fps()))

    def resizeGL(self, width, height):
        self.w, self.h = width, height
        self.paint_manager.resizeGL(width, height)

    def sizeHint(self):
        return QtCore.QSize(self.w, self.h)

    # Event methods
    # -------------
    def event(self, e):
        r = super(GalryWidget, self).event(e)
        if e.type() == QtCore.QEvent.Gesture:
            e.accept()
            gesture = e.gesture(QtCore.Qt.PinchGesture)
            self.pinchEvent(gesture)
            if gesture.state() == Qt.GestureStarted:
                self.mouse_blocked = True
            elif gesture.state() == Qt.GestureFinished:
                self.mouse_blocked = False
            return False
        return r

    def pinchEvent(self, e):
        self.user_action_generator.pinchEvent(e)
        self.process_interaction()

    def mousePressEvent(self, e):
        if self.mouse_blocked:
            return
        self.user_action_generator.mousePressEvent(e)
        self.process_interaction()

    def mouseReleaseEvent(self, e):
        if self.mouse_blocked:
            return
        self.user_action_generator.mouseReleaseEvent(e)
        self.process_interaction()

    def mouseDoubleClickEvent(self, e):
        if self.mouse_blocked:
            return
        self.user_action_generator.mouseDoubleClickEvent(e)
        self.process_interaction()

    def mouseMoveEvent(self, e):
        if self.mouse_blocked:
            return
        self.user_action_generator.mouseMoveEvent(e)
        self.process_interaction()

    def keyPressEvent(self, e):
        self.user_action_generator.keyPressEvent(e)
        self.process_interaction()
        # Close the application when pressing Q
        if e.key() == QtCore.Qt.Key_Q:
            if hasattr(self, 'window'):
                self.close_widget()

    def keyReleaseEvent(self, e):
        self.user_action_generator.keyReleaseEvent(e)

    def wheelEvent(self, e):
        self.user_action_generator.wheelEvent(e)
        self.process_interaction()

    def reset_action_generator(self):
        self.user_action_generator.reset()

    def leaveEvent (self, e):
        self.process_interaction(None)

    # Normalization methods
    # ---------------------
    def normalize_position(self, x, y):
        """Window coordinates ==> world coordinates."""
        if not hasattr(self.paint_manager, 'renderer'):
            return (0, 0)
        vx, vy = self.paint_manager.renderer.viewport
        x = -vx + 2 * vx * x / float(self.w)
        y = -(-vy + 2 * vy * y / float(self.h))
        return x, y

    def normalize_diff_position(self, x, y):
        """Normalize the coordinates of a difference vector between two
        points.
        """
        if not hasattr(self.paint_manager, 'renderer'):
            return (0, 0)
        vx, vy = self.paint_manager.renderer.viewport
        x = 2 * vx * x/float(self.w)
        y = -2 * vy * y/float(self.h)
        return x, y

    def normalize_action_parameters(self, parameters):
        """Normalize points in the action parameters object in the window
        coordinate system.

        Arguments:
          * parameters: the action parameters object, containing all
            variables related to user actions.

        Returns:
           * parameters: the updated parameters object with normalized
             coordinates.

        """
        parameters["mouse_position"] = self.normalize_position(\
                                                *parameters["mouse_position"])
        parameters["mouse_position_diff"] = self.normalize_diff_position(\
                                            *parameters["mouse_position_diff"])
        parameters["mouse_press_position"] = self.normalize_position(\
                                            *parameters["mouse_press_position"])
        parameters["pinch_position"] = self.normalize_position(\
                                            *parameters["pinch_position"])
        parameters["pinch_start_position"] = self.normalize_position(\
                                            *parameters["pinch_start_position"])
        return parameters

    # Signal methods
    # --------------
    def connect_events(self, arg1, arg2):
        """Makes a connection between a QT signal and an interaction event.

        The signal parameters must correspond to the event parameters.

        Arguments:
          * arg1: a QT bound signal or an interaction event.
          * arg2: an interaction event or a QT bound signal.

        """
        if type(arg1) == str:
            self.connect_event_to_signal(arg1, arg2)
        elif type(arg2) == str:
            self.connect_signal_to_event(arg1, arg2)

    def connect_signal_to_event(self, signal, event):
        """Connect a QT signal to an interaction event.

        The event parameters correspond to the signal parameters.

        Arguments:
          * signal: a QT signal.
          * event: an InteractionEvent string.

        """
        if signal is None:
            raise Exception("The signal %s is not defined" % signal)
        slot = lambda *args, **kwargs: \
                self.process_interaction(event, args, **kwargs)
        signal.connect(slot)

    def connect_event_to_signal(self, event, signal):
        """Connect an interaction event to a QT signal.

        The event parameters correspond to the signal parameters.

        Arguments:
          * event: an InteractionEvent string.
          * signal: a QT signal.

        """
        self.events_to_signals[event] = signal

    # Binding mode methods
    # --------------------
    def switch_interaction_mode(self):
        """Switch the interaction mode."""
        binding = self.binding_manager.switch()
        # set base cursor
        # self.interaction_manager.base_cursor = binding.base_cursor
        return binding

    def set_interaction_mode(self, mode):
        """Set the interaction mode.

        Arguments:
          * mode: either a class deriving from `Bindings` and which has been
            specified in `set_bindings`, or directly a `Bindings` instance.

        """
        binding = self.binding_manager.set(mode)
        # set base cursor
        # self.interaction_manager.base_cursor = binding.base_cursor
        return binding

    # Interaction methods
    # -------------------
    def get_current_action(self):
        """Return the current user action with the action parameters."""
        # get current action
        action = self.user_action_generator.action

        # get current key if the action was KeyPress
        key = self.user_action_generator.key

        # get key modifier
        key_modifier = self.user_action_generator.key_modifier

        # retrieve action parameters and normalize using the window size
        parameters = self.normalize_action_parameters(
                        self.user_action_generator.get_action_parameters())
        return action, key, key_modifier, parameters

    def get_current_event(self):
        """Return the current interaction event corresponding to the current
        user action."""
        # get the current interaction mode
        binding = self.binding_manager.get()

        # get current user action
        action, key, key_modifier, parameters = self.get_current_action()

        # get the associated interaction event
        event, param_getter = binding.get(action, key=key,
                                                  key_modifier=key_modifier)

        # get the parameter object by calling the param getter
        if param_getter is not None and parameters is not None:
            args = param_getter(parameters)
        else:
            args = None

        return event, args

    def set_current_cursor(self):
        cursor = self.interaction_manager.get_cursor()
        # if no cursor set, then use the default one in the current binding
        # mode
        if cursor is None:
            cursor = self.binding_manager.get().get_base_cursor()
        qcursor = get_cursor(cursor)
        if qcursor:
            self.setCursor(qcursor)

    def process_interaction(self, event=None, args=None, do_update=None):
        """Process user interaction.

        This method is called after each user action (mouse, keyboard...).
        It finds the right action associated to the command, then the event
        associated to that action.

        Arguments:
          * event=None: if None, the current event associated to the current
            user action is retrieved. Otherwise, an event can be directly
            passed here to force the trigger of any interaction event.
          * args=None: the arguments of the event if event is not None.

        """
        if event is None:
            # get current event from current user action
            event, args = self.get_current_event()

        if event == 'Animate' and self.block_refresh:
            return

        prev_event = self.interaction_manager.prev_event

        # handle interaction mode change
        if event == 'SwitchInteractionMode':
            binding = self.switch_interaction_mode()
            log_debug("Switching interaction mode to %s." % \
                binding.__class__.__name__)

        # process the interaction event
        self.interaction_manager.process_event(event, args)

        # raise a signal if there is one associated to the current event
        if event in self.events_to_signals:
            self.events_to_signals[event].emit(*args)

        # set cursor
        self.set_current_cursor()

        # clean current action (unique usage)
        self.user_action_generator.clean_action()

        # update the OpenGL view
        if do_update is None:
            do_update = (
                # (not isinstance(self, GalryTimerWidget)) and
                (event is not None or prev_event is not None))

        if do_update:
            self.updateGL()

    # Miscellaneous
    # -------------
    def save_image(self, file=None, update=True):
        """Save a screenshot of the widget in the specified file."""
        if file is None:
            file = "image.png"
        if update:
            self.updateGL()
        image = self.grabFrameBuffer()
        image.save(file, "PNG")

    def toggle_fullscreen(self):
        self.is_fullscreen = not self.is_fullscreen
        if self.is_fullscreen:
            if hasattr(self.window, 'showFullScreen'):
                self.window.showFullScreen()
        else:
            if hasattr(self.window, 'showNormal'):
                self.window.showNormal()

    def close_widget(self):
        self.user_action_generator.close()
        if hasattr(self, 'window'):
            if hasattr(self.window, 'close'):
                self.window.close()

    # Focus methods
    # -------------
    def focusOutEvent(self, event):
        self.user_action_generator.focusOutEvent(event)

class GalryTimerWidget(GalryWidget):
    timer = None

    """Special widget with periodic timer used to update the scene at
    regular intervals."""
    def initialize_timer(self, dt=1.):
        """Initialize the timer.

        This method *must* be called in the `initialize` method of the widget.

        Arguments:
          * dt=1.: the timer interval in seconds.

        """
        self.t = 0.
        self.dt = dt
        # start simulation after initialization completes
        self.timer = QtCore.QTimer()
        self.timer.setInterval(dt * 1000)
        self.timer.timeout.connect(self.update_callback)
        self.paint_manager.t = self.t

    def update_callback(self):
        """Callback function for the timer.

        Calls `paint_manager.update_callback`, so this latter method should be
        implemented in the paint manager. The attribute `self.t` is
        available here and in the paint manager.

        """
        self.t = timeit.default_timer() - self.t0
        self.process_interaction('Animate', (self.t,))

    def start_timer(self):
        """Start the timer."""
        if self.timer:
            self.t0 = timeit.default_timer()
            self.timer.start()

    def stop_timer(self):
        """Stop the timer."""
        if self.timer:
            self.timer.stop()

    def showEvent(self, e):
        """Called when the window is shown (for the first time or after
        minimization). It starts the timer."""
        # start simulation when showing window
        self.start_timer()

    def hideEvent(self, e):
        """Called when the window is hidden (e.g. minimized). It stops the
        timer."""
        # stop simulation when hiding window
        self.stop_timer()

# Basic widgets helper functions and classes
# ------------------------------------------
def create_custom_widget(bindings=None,
                         antialiasing=False,
                         constrain_ratio=False,
                         constrain_navigation=False,
                         activate_help=True,
                         activate_grid=False,
                         show_grid=False,
                         display_fps=False,
                         activate3D=False,
                         animation_interval=None,
                         momentum=False,
                         autosave=None,
                         getfocus=True,
                         figure=None,
                        **companion_classes):
    """Helper function to create a custom widget class from various parameters.

    Arguments:
      * bindings=None: the bindings class, instance, or a list of those.
      * antialiasing=False: whether to activate antialiasing or not. It can
        hurt performance.
      * constrain_ratio=False: if True, the ratio is 1:1 at all times.
      * constrain_navigation=True: if True, the viewbox cannot be greater
        than [-1,1]^2 by default (but it can be customized in
        interactionmanager.MAX_VIEWBOX).
      * display_fps=False: whether to display the FPS.
      * animation_interval=None: if not None, a special widget with automatic
        timer update is created. This variable then refers to the time interval
        between two successive updates (in seconds).
      * **companion_classes: keyword arguments with the companion classes.

    """
    if momentum and animation_interval is None:
        animation_interval = .01

    # use the GalryTimerWidget if animation_interval is not None
    if animation_interval is not None:
        baseclass = GalryTimerWidget
    else:
        baseclass = GalryWidget

    if bindings is None:
        bindings = []
    if type(bindings) != list and type(bindings) != tuple:
        bindings = [bindings]

    # create the custom widget class
    class MyWidget(baseclass):
        """Automatically-created Galry widget."""
        def __init__(self):
            # antialiasing
            if QGLFormat is not None:
                format = QGLFormat()
            else:
                format = None
            if antialiasing:
                if hasattr(format, 'setSampleBuffers'):
                    format.setSampleBuffers(True)
            super(MyWidget, self).__init__(format=format, autosave=autosave,
                getfocus=getfocus)

        def initialize(self):
            if figure:
                figure.widget = self
            self.set_bindings(*bindings)
            self.set_companion_classes(**companion_classes)
            self.constrain_ratio = constrain_ratio
            self.constrain_navigation = constrain_navigation
            self.activate_help = activate_help
            self.activate_grid = activate_grid
            self.show_grid = show_grid
            self.activate3D = activate3D
            self.momentum = momentum
            self.display_fps = display_fps
            self.initialize_companion_classes()
            if animation_interval is not None:
                self.initialize_timer(dt=animation_interval)

    return MyWidget

class AutodestructibleWindow(QtGui.QMainWindow):
    """Special QT window that can be destroyed automatically after a given
    timeout. Useful for automatic debugging or benchmarking."""
    autodestruct = None

    def __init__(self, **kwargs):
        super(AutodestructibleWindow, self).__init__()
        # This is important in interaction sessions: it allows the widget
        # to clean everything up as soon as we close the window (otherwise
        # it is just hidden).
        self.setAttribute(QtCore.Qt.WA_DeleteOnClose)
        self.initialize(**kwargs)

    def set_autodestruct(self, autodestruct=None):
        # by default, use global variable
        if autodestruct is None:
            # use the autodestruct option in command line by default
            autodestruct = "autodestruct" in sys.argv
            if autodestruct is False:
                global AUTODESTRUCT
                autodestruct = AUTODESTRUCT
        # option for autodestructing the window after a fixed number of
        # seconds: useful for automatic testing
        if autodestruct is True:
            # 3 seconds by default, if True
            global DEFAULT_AUTODESTRUCT
            autodestruct = DEFAULT_AUTODESTRUCT
        if autodestruct:
            log_info("window autodestruction in %d second(s)" % (autodestruct / 1000.))
        self.autodestruct = autodestruct

    def initialize(self, **kwargs):
        pass

    def kill(self):
        if self.autodestruct:
            self.timer.stop()
            self.close()

    def showEvent(self, e):
        if self.autodestruct:
            self.timer = QtCore.QTimer()
            self.timer.setInterval(self.autodestruct)
            self.timer.setSingleShot(True)
            self.timer.timeout.connect(self.kill)
            self.timer.start()

def create_basic_window(widget=None, size=None, position=(20, 20),
                        autodestruct=None,
                        toolbar=False):
    """Create a basic QT window with a Galry widget inside.

    Arguments:
      * widget: a class or instance deriving from GalryWidget.
      * size=None: the size of the window as a tuple (width, height).
      * position=(100, 100): the initial position of the window on the screen,
        in pixels (x, y).
      * autodestruct=None: if not None, it is the time, in seconds, before the
        window closes itself.

    """
    class BasicWindow(AutodestructibleWindow):
        """Automatically-created QT window."""
        def initialize(self, widget=widget, size=size, position=position,
                       autodestruct=autodestruct):
            """Create a basic window to display a single widget.

            Arguments:
              * widget: a class or instance deriving from GalryWidget.
              * size=None: the size of the window as a tuple (width, height).
              * position=(100, 100): the initial position of the window on the screen,
                in pixels (x, y).
              * autodestruct=None: if not None, it is the time, in seconds, before the
                window closes itself.

            """
            self.set_autodestruct(autodestruct)
            # default widget
            if widget is None:
                widget = GalryWidget()
            # if widget is not an instance of GalryWidget, maybe it's a class,
            # then try to instanciate it
            if not isinstance(widget, GalryWidget):
                widget = widget()
            widget.window = self
            # create widget
            self.widget = widget
            if toolbar:
                self.add_toolbar()
            # if size is None:
                # size = self.widget.w, self.widget.h
            if size is not None:
                self.widget.w, self.widget.h = size
            self.setCentralWidget(self.widget)
            self.setWindowTitle("Galry")
            self.move(*position)
            # ensure the main window size is adjusted so that the widget size
            # is equal to the specified size
            self.resize(self.sizeHint())
            self.show()

        def toggle_toolbar(self):
            self.toolbar.setVisible(not self.toolbar.isVisible())#not )

        def add_toolbar(self):
            """Add navigation toolbar"""

            # reset
            reset_action = QtGui.QAction("Reset view (R)", self)
            reset_action.setIcon(get_icon('home'))
            self.widget.connect_events(reset_action.triggered, 'Reset')

            # show grid
            grid_action = QtGui.QAction("Show grid (G)", self)
            grid_action.setIcon(get_icon('grid'))
            self.widget.connect_events(grid_action.triggered, 'Grid')

            # fullscreen
            fullscreen_action = QtGui.QAction("Fullscreen (F)", self)
            fullscreen_action.setIcon(get_icon('fullscreen'))
            self.widget.connect_events(fullscreen_action.triggered, 'Fullscreen')

            # save image
            save_action = QtGui.QAction("Save image (S)", self)
            save_action.setIcon(get_icon('save'))
            save_action.setShortcut("S")
            save_action.triggered.connect(self.save)

            toolbar_action = QtGui.QAction("Toggle toolbar visibility (T)", self)
            toolbar_action.setIcon(get_icon('toolbar'))
            toolbar_action.setShortcut("T")
            toolbar_action.triggered.connect(self.toggle_toolbar)
            # self.toolbar_action = toolbar_action

            # help
            help_action = QtGui.QAction("Show help (H)", self)
            help_action.setIcon(get_icon('help'))
            self.widget.connect_events(help_action.triggered, 'Help')

            # exit
            exit_action = QtGui.QAction("Exit (Q)", self)
            exit_action.setIcon(get_icon('exit'))
            exit_action.triggered.connect(self.close)

            # add toolbar
            mytoolbar = QtGui.QToolBar(self.widget)
            mytoolbar.setIconSize(QtCore.QSize(32, 32))

            for action in [reset_action, grid_action, fullscreen_action,
                toolbar_action, save_action, help_action, exit_action]:
                self.addAction(action)
                mytoolbar.addAction(action)

            mytoolbar.setStyleSheet("""
            QToolBar, QToolButton
            {
                background: #000000;
                border-color: #000000;
                color: #ffffff;
            }
            QToolButton
            {
                margin-left: 5px;
            }
            QToolButton:hover
            {
                background: #2a2a2a;
            }
            """)
            mytoolbar.setMovable(False)
            mytoolbar.setFloatable(False)
            self.toolbar = mytoolbar
            self.addToolBar(mytoolbar)

        def save(self):
            """Open a file dialog and save the current image in the specified
            PNG file."""
            initial_filename = 'screen'
            existing = filter(lambda f: f.startswith(initial_filename), os.listdir('.'))
            i = 0
            if existing:
                for f in existing:
                    r = re.match('screen([0-9]*).png', f)
                    i = max(i, int(r.groups()[0]))
                i += 1
                # if last:
                    # last = int(last)
                    # i = last + 1
            filename, _ = QtGui.QFileDialog.getSaveFileName(self,
                "Save the current view in a PNG image",
                initial_filename + str(i) + '.png',
                # '*.png',
                # '*.png',
                # QtGui.QFileDialog.AnyFile,
                )
            if filename:
                self.widget.save_image(str(filename))

        def closeEvent(self, e):
            """Clean up memory upon closing."""
            self.widget.user_action_generator.close()
            self.widget.paint_manager.cleanup()
            super(BasicWindow, self).closeEvent(e)

        def contextMenuEvent(self, e):
            return

    return BasicWindow

def show_basic_window(widget_class=None, window_class=None, size=None,
            position=(20, 20), autodestruct=None, toolbar=False, **kwargs):
    """Create a custom widget and/or window and show it immediately.

    Arguments:
      * widget_class=None: the class deriving from GalryWidget.
      * window_class=None: the window class, deriving from `QMainWindow`.
      * size=None: the size of the window as a tuple (width, height).
      * position=(100, 100): the initial position of the window on the screen,
        in pixels (x, y).
      * autodestruct=None: if not None, it is the time, in seconds, before the
        window closes itself.
      * **kwargs: keyword arguments with the companion classes and other
        parameters that are passed to `create_custom_widget`.

    """
    # default widget class
    if widget_class is None:
        widget_class = create_custom_widget(**kwargs)
    # defaut window class
    if window_class is None:
        window_class = create_basic_window(widget_class, size=size,
            position=position, autodestruct=autodestruct, toolbar=toolbar,
            )
    # create and show window
    return show_window(window_class)

max_stars_repo_name: fiath/test

max_issues_repo_name: fiath/test

max_forks_repo_name: fiath/test


text:
import argparse
import numpy as np
from numpy import *
import matplotlib.pyplot as plt
import networkx as nx
import grnular.data.gen_data as gen_data
from grnular.utils.metrics import report_metrics
import sys, copy, pickle
from sklearn.covariance import GraphicalLassoCV
from sklearn.covariance import GraphicalLasso
import sklearn
TRAIN=True

parser = argparse.ArgumentParser(description='Classification of different cell types as well as recovering the gene regulatory network of RNA seq data: using SERGIO simulator for training')

#****************** general parameters
parser.add_argument('--K_train', type=int, default=5, #2, #1000,
                    help='Num of training examples for a fixed D')
parser.add_argument('--K_valid', type=int, default=5, #5, #2, #1000,
                    help='Number of valid examples for a fixed D ')
parser.add_argument('--K_test', type=int, default=100, #100, #10,
                    help='Number of testing examples for a fixed D')
parser.add_argument('--D', type=int, default=100, #1000,
                    help='Number of genes ')
parser.add_argument('--C', type=int, default=9,
                    help='different cell types, target variable')
parser.add_argument('--sparsity', type=float, default=0.3, #0.2,
                    help='sparsity of erdos-renyi graph')
parser.add_argument('--DATA_METHOD', type=str,  default='sim_expt', #
                    help='expt details in draft: random/syn_same_precision, sim_expt1=DS1 and sim_expt2=Custom, sim_expt3=GRN')
parser.add_argument('--DATA_TYPE', type=str,  default='clean', #'case2',
                    help='expt details in draft: sim_exp1 clean/noisy')
parser.add_argument('--alpha_l1', type=float, default=0.1, #0.2,
                    help='L1 penalty for graphical lasso')
parser.add_argument('--mode', type=str,  default='cd', #'lars'
                    help='glasso mode: cd or lars or cv')

# SERGIO simulator parameters
parser.add_argument('--DATA_NAME', type=str,  default='CUSTOM', #'DS1',
                    help='expt details in draft: DS1, DS2, DS3, CUSTOM')
parser.add_argument('--POINTS_PER_CLASS', type=int, default=2000,# NOTE: try 2000
                    help='cells per class type')
parser.add_argument('--SAMPLING_STATE', type=int, default=1, #1,
                    help='num of simulations')
parser.add_argument('--NOISE_PARAMS', type=float, default=0.1, #1,
                    help='num of noise params')
parser.add_argument('--DECAYS', type=float, default=5, #0.8,#0.8,
                    help='decay params')
parser.add_argument('--NOISE_TYPE', type=str,  default='dpd', #'dpd',
                    help='different noise types: "dpd", “sp”, “spd”')
parser.add_argument('--SHARED_COOP_STATE', type=int, default=2, #1,
                    help='shared coop state')
parser.add_argument('--pcr_low_min', type=float, default=0.2, #1,
                    help='production cell rate: low expression range')
parser.add_argument('--pcr_low_max', type=float, default=0.5, #1,
                    help='production cell rate: low expression range')
parser.add_argument('--pcr_high_min', type=float, default=0.7, #1,
                    help='production cell rate: high expression range')
parser.add_argument('--pcr_high_max', type=float, default=1, #1,
                    help='production cell rate: high expression range')
parser.add_argument('--Kij_min', type=float, default=1, #1,
                    help='Interaction strengths Kij min')
parser.add_argument('--Kij_max', type=float, default=5, #1,
                    help='Interaction strengths Kij max')
parser.add_argument('--ratio_MR', type=float, default=0.1, #1,
                    help='number of master regulators ~ ratio_MR * D')
parser.add_argument('--connect_TF_prob', type=float, default=0.2, #1,
                    help='probability of connecting master regulators')
parser.add_argument('--USE_TF_NAMES', type=str,  default='no',# 'yes'
                    help='use transcription factors to improve prediction: will use in general')
# SERGIO technical noise parameters
parser.add_argument('--ADD_TECHNICAL_NOISE', type=str,  default='yes',# 'no'
                    help='add technical noise on the saved clean data')
parser.add_argument('--dropout_shape', type=float, default=6.5, #1,
                    help='SERGIO dropout param: shape, higher -> less dropout')
parser.add_argument('--dropout_percentile', type=float, default=82, #1,
                    help='SERGIO dropout param: percentile, lower -> less dropout')

args = parser.parse_args()
print(args)
print('\n')

def get_args_str(dict1):
    args_str = ''
    for i, (k, v) in enumerate(dict1.items()):
#        print(k , item)
        if k in ['C', 'GLAD_LOSS', 'MODEL_SELECT', 'SUB_METHOD']:
            args_str = args_str+str(k)+str(v)+'_'
    return args_str

args_str = get_args_str(vars(args))

def fit_graphical_lasso(data, PREDICT_TF=False, BEELINE=False):
    print('FITTING a graphical lasso')
    # #############################################################################
    res = []
    typeS = 'mean'
    print('Using ', typeS, ' scaling')
    for i, d in enumerate(data):
        X, y, theta_true, master_regulators = d
        Xc = normalizing_data(X, typeS)
        print('\nGRAPHICAL Lasso: TRAIN data batch : ', i, ' total points = ', X.shape[0])
        if args.USE_TF_NAMES=='yes' and PREDICT_TF:
            res.append(helper_graphical_lasso(Xc, theta_true, tf_names=master_regulators))
        else:
            res.append(helper_graphical_lasso(Xc, theta_true))
    res_mean = np.mean(np.array(res), 0)
    res_std  = np.std(np.array(res).astype(np.float64), 0)
    res_mean = ["%.3f" %x for x in res_mean]
    res_std  = ["%.3f" %x for x in res_std]
    res_dict = {} # name: [mean, std]
    for i, _name in enumerate(['FDR', 'TPR', 'FPR', 'SHD', 'nnz_true', 'nnz_pred', 'precision', 'recall', 'Fb', 'aupr', 'auc']): # dictionary
        res_dict[_name]= [res_mean[i], res_std[i]]#mean std
    if PREDICT_TF:
        print('\nAvg GLASSO-TF: FDR, ,TPR, ,FPR, ,SHD, ,nnz_true, ,nnz_pred, ,precision, ,recall, ,Fb, ,aupr, ,auc, ')
    else:
        print('\nAvg GLASSO: FDR, ,TPR, ,FPR, ,SHD, ,nnz_true, ,nnz_pred, ,precision, ,recall, ,Fb, ,aupr, ,auc, ')
    mean_std = [[rm, rs] for rm, rs in zip(res_mean, res_std)]
    flat_list = [item for ms in mean_std for item in ms]
    print('%s' % ', '.join(map(str, flat_list)))
    #print(*sum(list(map(list, zip(res_mean, res_std))), []), sep=', ')
    return

def normalizing_data(X, typeS='log'):
    if typeS == 'mean':
        #print('Centering and scaling the input data...')
        scaledX = X - X.mean(axis=0)
        scaledX = scaledX/X.std(axis=0)
        # NOTE: replacing all nan's by 0, as sometimes in dropout the complete column
        # goes to zero
        scaledX = convert_nans_to_zeros(scaledX)
    elif typeS == 'log':
        scaledX = np.log(X+1)
    else:
        print('Check the valid scaling')
    return scaledX

def convert_nans_to_zeros(X):
    where_are_nans = isnan(X)
    X[where_are_nans] = 0
    return X

def get_PSD_matrix(A, u=1):
    smallest_eigval = np.min(np.linalg.eigvals(A))
    # Just in case : to avoid numerical error in case a epsilon complex component present
    smallest_eigval = smallest_eigval.real
    # making the min eigenvalue as 1
    target_precision_mat = A + np.eye(A.shape[-1])*(u - smallest_eigval)
    #print('CHEKKK: smallest eigen? = ', np.min(np.linalg.eigvals(target_precision_mat)))
    return target_precision_mat

def postprocess_tf(prec, tf_names):
#    print('Postprocesing for TF NAMES')
    # remove all the edges whose at least one of the vertices is not in tf_names
    # zeroing the diagonal to get adj matrix
    tf_names = ['G'+str(n) for n in tf_names]
    np.fill_diagonal(prec, 0)
    G_pred = nx.from_numpy_matrix(prec)
    #mapping = {n: 'G'+str(n) for n in range(args.D)}
    mapping = {n: 'G'+str(n) for n in range(prec.shape[-1])}
    G_pred = nx.relabel_nodes(G_pred, mapping)
    set_tf_names = tf_names
    edges = copy.deepcopy(G_pred.edges)
    for e in edges:
        if not (set(e) & set(tf_names)):
            # remove the edge
            G_pred.remove_edge(*e[:2])
    prec = get_PSD_matrix(nx.adj_matrix(G_pred).todense())
    return np.array(prec)

def helper_graphical_lasso(X, theta_true, tf_names=[]):
    # Estimate the covariance
    if args.mode == 'cv':
        model = GraphicalLassoCV()
    else:
        model = GraphicalLasso(alpha=args.alpha_l1, mode=args.mode, tol=1e-7, enet_tol=1e-6,
                                max_iter=100, verbose=False, assume_centered=False)
    model.fit(X)
#    cov_ = model.covariance_
    prec_ = model.precision_
    if args.USE_TF_NAMES=='yes' and len(tf_names) != 0:
        prec_ = postprocess_tf(prec_, tf_names)
    recovery_metrics = report_metrics(np.array(theta_true), prec_)
    print('GLASSO: FDR, TPR, FPR, SHD, nnz_true, nnz_pred, precision, recall, Fb, aupr, auc')
    print('GLASSO: TEST: Recovery of true theta: ', *np.around(recovery_metrics, 3))
    return list(recovery_metrics)

def get_filepath():
    dict1 = vars(args)
    filename = ''
    abbrv_dict = {'K_train': 'KTr', 'K_valid': 'KVa', 'K_test': 'KTe', 'D': 'D', 'C':'C',
                  'sparsity': 'Sp', 'DATA_TYPE':'Dt', 'POINTS_PER_CLASS': 'ppc',
                  'SAMPLING_STATE': 'SS', 'NOISE_PARAMS': 'NP', 'DECAYS': 'De',
                  'NOISE_TYPE': 'NT', 'SHARED_COOP_STATE': 'SCS', 'pcr_low_min': 'pcrln',
                  'pcr_low_max': 'pcrlx', 'pcr_high_min': 'pcrhn', 'pcr_high_max': 'pcrhx',
                  'Kij_min': 'kmin', 'Kij_max': 'kmax', 'ratio_MR': 'rMR',
                  'connect_TF_prob': 'TFp'}
    for k in abbrv_dict.keys():
        v = dict1[k]
        filename = filename+str(abbrv_dict[k])+str(v)+'_'

    SAVEPATH = 'grnular/data/saved_data/'
    FILEPATH = SAVEPATH + filename + '.pickle'
    print('Filepath: ', FILEPATH)
    return FILEPATH

def load_saved_data():
    FILEPATH = get_filepath()
    with open(FILEPATH, 'rb') as handle:
        data = pickle.load(handle)
    return data

def main():
    print('Reading the input data: Single cell RNA: M(samples) x D(genes) & corresponding C(targets)')
    train_data, valid_data, test_data = load_saved_data()

    if args.ADD_TECHNICAL_NOISE == 'yes':
        print('adding technical noise')
#        train_data = gen_data.add_technical_noise(args, train_data)
#        valid_data = gen_data.add_technical_noise(args, valid_data)
        test_data = gen_data.add_technical_noise(args, test_data)

    # Fitting a graphical lasso
    fit_graphical_lasso(test_data)
    print('Using TF NAMES')
    fit_graphical_lasso(test_data, PREDICT_TF=True)
    print('\nExpt Done')
    return

if __name__=="__main__":
    main()

ext: py

reason: import numpy,from numpy,import networkx


text:
import unittest
import numpy as np

from flowutils import transforms

class TransformsTestCase(unittest.TestCase):
    def setUp(self):
        self.test_data_range = np.linspace(0.0, 1000.0, 10001)

    @staticmethod
    def test_logicle_range():
        """Test a range of input values"""
        data_in = np.array([-10, -5, -1, 0, 0.3, 1, 3, 10, 100, 1000])
        correct_output = np.array(
            [
                0.067574,
                0.147986,
                0.228752,
                0.25,
                0.256384,
                0.271248,
                0.312897,
                0.432426,
                0.739548,
                1.0
            ]
        )

        # noinspection PyProtectedMember
        data_out = transforms._logicle(data_in, t=1000, m=4.0, w=1.0, a=0)

        np.testing.assert_array_almost_equal(data_out, correct_output, decimal=6)

    def test_inverse_logicle_transform(self):
        xform_data = transforms.logicle(
            self.test_data_range.reshape(-1, 1),
            [0],
            t=10000,
            w=0.5,
            m=4.5,
            a=0
        )
        x = transforms.logicle_inverse(
            xform_data,
            [0],
            t=10000,
            w=0.5,
            m=4.5,
            a=0
        )

        np.testing.assert_array_almost_equal(self.test_data_range, x[:, 0], decimal=10)

    @staticmethod
    def test_asinh_range():
        """Test a range of input values"""
        data_in = np.array([-10.0, -5.0, -1.0, 0.0, 0.3, 1.0, 3.0, 10.0, 100.0, 1000.0], dtype=np.float)
        data_in = data_in.reshape((-1, 1))
        correct_output = np.array(
            [[
                -0.200009,
                -0.139829,
                -0.000856,
                0.2,
                0.303776,
                0.400856,
                0.495521,
                0.600009,
                0.8,
                1.0
            ]]
        ).reshape((-1, 1))

        data_out = transforms.asinh(data_in, channel_indices=0, t=1000, m=4.0, a=1.0)

        np.testing.assert_array_almost_equal(data_out, correct_output, decimal=6)

    def test_inverse_asinh_transform(self):
        xform_data = transforms.asinh(
            self.test_data_range.reshape(-1, 1),
            [0],
            t=10000,
            m=4.5,
            a=0
        )
        x = transforms.asinh_inverse(
            xform_data,
            [0],
            t=10000,
            m=4.5,
            a=0
        )

        np.testing.assert_array_almost_equal(self.test_data_range, x[:, 0], decimal=10)

    @staticmethod
    def test_hyperlog_range():
        """Test a range of input values"""
        data_in = np.array([-10, -5, -1, 0, 0.3, 1, 3, 10, 100, 1000])
        correct_output = np.array(
            [
                0.08355406,
                0.15586819,
                0.2294768,
                0.25,
                0.25623887,
                0.2705232,
                0.30909185,
                0.41644594,
                0.73187469,
                1.
            ]
        )

        # noinspection PyProtectedMember
        data_out = transforms._hyperlog(data_in, t=1000, m=4.0, w=1.0, a=0)

        np.testing.assert_array_almost_equal(data_out, correct_output, decimal=6)

    def test_inverse_hyperlog_transform(self):
        xform_data = transforms.hyperlog(
            self.test_data_range.reshape(-1, 1),
            [0],
            t=10000,
            w=0.5,
            m=4.5,
            a=0
        )
        x = transforms.hyperlog_inverse(
            xform_data,
            [0],
            t=10000,
            w=0.5,
            m=4.5,
            a=0
        )

        np.testing.assert_array_almost_equal(self.test_data_range, x[:, 0], decimal=10)

    @staticmethod
    def test_log_range():
        """Test a range of input values"""
        data_in = np.array(
            [-1., 0., 0.5, 1., 10., 100., 1000., 1023., 10000., 100000, 262144],
            dtype=np.float
        )
        data_in = data_in.reshape((-1, 1))
        correct_output = np.array(
            [[
                np.nan,
                -np.inf,
                0.139794,
                0.2,
                0.4,
                0.6,
                0.8,
                0.801975,
                1.0,
                1.2,
                1.283708
            ]]
        ).reshape((-1, 1))

        with np.errstate(divide='ignore', invalid='ignore'):
            data_out = transforms.log(data_in, channel_indices=0, t=10000, m=5.0)

        np.testing.assert_array_almost_equal(data_out, correct_output, decimal=6)

    def test_inverse_log_transform(self):
        with np.errstate(divide='ignore'):
            xform_data = transforms.log(
                self.test_data_range.reshape(-1, 1),
                [0],
                t=10000,
                m=4.5
            )
        x = transforms.log_inverse(
            xform_data,
            [0],
            t=10000,
            m=4.5
        )

        np.testing.assert_array_almost_equal(self.test_data_range, x[:, 0], decimal=10)

max_stars_repo_name: DavHau/FlowUtils

max_issues_repo_name: DavHau/FlowUtils

max_forks_repo_name: DavHau/FlowUtils

reason: import numpy


text:
# Copyright 2018 Amazon.com, Inc. or its affiliates. All Rights Reserved.
#
# Licensed under the Apache License, Version 2.0 (the "License").
# You may not use this file except in compliance with the License.
# A copy of the License is located at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# or in the "license" file accompanying this file. This file is distributed
# on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either
# express or implied. See the License for the specific language governing
# permissions and limitations under the License.

import mxnet as mx
import numpy as np

class MetricAttentiveScheduler(mx.lr_scheduler.LRScheduler):
    r"""
    This scheduler decreases the learning rate based on the value of some
    validation metric to be optimized (maximized or minimized). The value
    of such metric is provided by calling the `step` method on the scheduler.
    A `patience` parameter must be provided, and the scheduler will reduce
    the learning rate if no improvement in the metric is done before
    `patience` observations of the metric.

    Examples:

        `patience = 0`: learning rate will decrease at every call to
        `step`, regardless of the metric value

        `patience = 1`: learning rate is reduced as soon `step` is called
        with a metric value which does not improve over the best encountered

        `patience = 10`: learning rate is reduced if no improvement in the
        metric is recorded in 10 successive calls to `step`

    Parameters
    ----------
    objective
        String, can either be `"min"` or `"max"`
    patience
        The patience to observe before reducing the learning rate, nonnegative integer.
    base_lr
        Initial learning rate to be used.
    decay_factor
        Factor (between 0 and 1) by which to decrease the learning rate.
    min_lr
        Lower bound for the learning rate, learning rate will never go below `min_lr`
    """

    def __init__(
        self,
        objective: str,
        patience: int,
        base_lr: float = 0.01,
        decay_factor: float = 0.5,
        min_lr: float = 0.0,
    ) -> None:

        assert base_lr > 0, f"base_lr should be positive, got {base_lr}"

        assert (
            base_lr > min_lr
        ), f"base_lr should greater than min_lr, {base_lr} <= {min_lr}"

        assert (
            0 < decay_factor < 1
        ), f"decay_factor factor should be between 0 and 1, got {decay_factor}"

        assert patience >= 0, f"patience should be nonnegative, got {patience}"

        assert objective in [
            "min",
            "max",
        ], f"objective should be 'min' or 'max', got {objective}"

        super(MetricAttentiveScheduler, self).__init__(base_lr=base_lr)

        self.decay_factor = decay_factor
        self.patience = patience
        self.objective = objective
        self.min_lr = min_lr
        self.best_metric = np.Inf if objective == "min" else -np.Inf
        self.prev_change = 0
        self.epoch_no = 0
        self.curr_lr = None

    def __call__(self, num_update: int) -> float:
        if self.curr_lr is None:
            self.curr_lr = self.base_lr
        assert self.curr_lr is not None

        return self.curr_lr

    def step(self, metric_value: float) -> bool:
        """
        Inform the scheduler of the new value of the metric that is being
        optimized. This method should be invoked at regular intervals (e.g.
        at the end of every epoch, after computing a validation score).

        Parameters
        ----------
        metric_value
            Value of the metric that is being optimized.

        Returns
        -------
        bool value indicating, whether to continue training
        """
        if self.curr_lr is None:
            self.curr_lr = self.base_lr
        assert self.curr_lr is not None

        metric_improved = (
            self.objective == "min" and metric_value < self.best_metric
        ) or (self.objective == "max" and metric_value > self.best_metric)

        if metric_improved:
            self.best_metric = metric_value
            self.prev_change = self.epoch_no

        if (
            self.epoch_no - self.prev_change >= self.patience
            or not np.isfinite(metric_value)
        ):
            if self.curr_lr == self.min_lr:
                return False
            self.curr_lr = max(self.min_lr, self.decay_factor * self.curr_lr)
            self.prev_change = self.epoch_no

        self.epoch_no += 1
        return True


text:
#!/usr/bin/env python3
import numpy as np

def fix_pixels(image, bad_pixels):
    bad_indices = np.nonzero(bad_pixels)
    mean = np.zeros(bad_indices[0].size, dtype=np.float32)
    mean_good = np.zeros(bad_indices[0].size, dtype=np.float32)
    norm = np.zeros(bad_indices[0].size, dtype=np.int32)
    shifts = ((-1,-1),
              (-1, 0),
              (-1, 1),
              ( 0, 1),
              ( 1, 1),
              ( 1, 0),
              ( 1,-1),
              ( 0,-1))

    for di, dj in shifts:
        I = np.clip(bad_indices[0] + di, 0, image.shape[0]-1)
        J = np.clip(bad_indices[1] + dj, 0, image.shape[1]-1)
        px = image[I, J]
        bmask = bad_pixels[I, J]
        mean += px
        mean_good += np.where(bmask, 0, px)
        norm += (~bmask).astype(np.int32)
    mean /= len(shifts)
    mean_good /= np.maximum(1, norm)
    image[bad_indices[0], bad_indices[1]] = np.where(norm <= len(shifts)/2, mean, mean_good)

if __name__ == "__main__":
    import matplotlib.pyplot as plt
    from medpy.filter.smoothing import anisotropic_diffusion
    from scipy.ndimage.filters import gaussian_filter
    from map_envi import load

    #imt = np.load('tnsr.npy')[3000:4000,4000:5000,0]
    #bdt = np.load('bad_data.npy')[3000:4000,4000:5000]
    #print(imt.min(), imt.max(), imt.mean(), imt.std())

    # im = np.zeros((1000, 1000), dtype=np.float32)
    # bd = np.zeros((1000, 1000), dtype=np.bool)
    # for fn in ('Sigma0_IW1_VH_mst_22Jun2018',
    #            'Sigma0_IW1_VH_slv1_04Jul2018',
    #            'Sigma0_IW1_VV_mst_22Jun2018',
    #            'Sigma0_IW1_VV_slv2_04Jul2018'):
    #     patch, _ = load(fn)
    #     patch = patch[10000:11000,11000:12000].astype(np.float32)
    #     bd |= (patch <= 1e-6) | (patch > 10)
    #     im += patch

    # im *= 0.25
    #im = np.where(bd, 0, np.log10(np.maximum(1e-6, im)) )

    im1 = load('Sigma0_IW1_VV_mst_22Jun2018')[0][7000:13000,7000:13000].astype(np.float32)
    im2 = load('Sigma0_IW1_VH_mst_22Jun2018')[0][7000:13000,7000:13000].astype(np.float32)

    im = np.hypot(im1, im2)
    bd = (im <= 1e-6) | (im > 10)

    im = np.arcsin(im2 / im)

    im = np.maximum(1e-6, im)
    imor = im.copy()

    fix_pixels(im, bd)

    u = gaussian_filter(im, 4)
    u = anisotropic_diffusion(im, 50, 20, 0.25, option=1)
    #u = im
    G1 = np.hypot(*np.gradient(im))
    G2 = np.hypot(*np.gradient(u))

    plt.subplot(2, 2, 1), plt.imshow(imor, cmap='gray')
    plt.xticks([]), plt.yticks([])

    plt.subplot(2, 2, 2), plt.imshow(u , cmap='gray')
    plt.xticks([]), plt.yticks([])

    plt.subplot(2, 2, 3), plt.imshow(G1, cmap='gray')
    plt.xticks([]), plt.yticks([])

    plt.subplot(2, 2, 4), plt.imshow(G2, cmap='gray')
    plt.xticks([]), plt.yticks([])

    plt.subplots_adjust(bottom=0., right=1, top=1)

    plt.show()

max_stars_repo_name: geobeyond/AI-GEO

max_issues_repo_name: geobeyond/AI-GEO

max_forks_repo_name: geobeyond/AI-GEO

reason: import numpy,from scipy


text:
import numpy as np

import theano
import theano.tensor as T
import lasagne

from theano.ifelse import ifelse

from neuralnet import iterate_minibatches

from collections import OrderedDict
import time
import random

EXTRA_INFO = True
DEFAULT_ADAPTIVE = False

# Set up parameters for Streaming SVRG
SMOOTHNESS = 1

class StreamingSVRGOptimizer:
    def __init__(self, m, learning_rate, adaptive=DEFAULT_ADAPTIVE, non_uniform_prob=True):
        self.m = m
        self.learning_rate = learning_rate
        # Adaptive is if we use line search to dynamically decide the learning rate.
        self.adaptive = adaptive
        self.non_uniform_prob = non_uniform_prob
        self.counted_gradient = theano.shared(0)
        self.smoothness = SMOOTHNESS

    def minimize(self, loss, params, X_train, Y_train, X_test, y_test, input_var, target_var, X_val, Y_val, n_epochs=1000, batch_size=100, output_layer=None ):
        self.input_var = input_var
        self.target_var = target_var

        flog = open("data/log.txt",'w')

        num_batches = X_train.shape[0] / batch_size
        n = num_batches

        if EXTRA_INFO:
            flog.write("Learning Rate:{:.2f}\n".format(self.learning_rate))
            flog.write("Adaptive:{:.2f}\n".format(self.adaptive))
            flog.write("Non Uniform Prob of mini batch:{:.2f}\n".format(self.non_uniform_prob))

        self.L = theano.shared(np.cast['float32'](1. / self.learning_rate))
        #self.Ls = [theano.shared(np.cast['float32'](1. / self.learning_rate)) for _  in range(num_batches)]
        self.Ls = [1. / self.learning_rate for _  in range(num_batches)]

        # params are all trainable parameters
        w_updates, mu_updates = self.make_updates(loss, params)

        train_mu = theano.function([self.input_var, self.target_var], loss, updates=mu_updates)
        train_w = theano.function([self.input_var, self.target_var], loss, updates=w_updates)

        prediction = lasagne.layers.get_output(output_layer, deterministic=True)
        acc_fn =      T.mean(T.eq(T.argmax(prediction, axis=1), target_var),dtype=theano.config.floatX)
        train_acc_fn = theano.function([self.input_var, self.target_var], acc_fn)
        val_fn = theano.function([self.input_var, self.target_var], [loss, acc_fn])

        train_error = []
        validation_error = []
        acc_train = []
        acc_val = []
        test_error = []
        acc_test = []
        times = []

        print "NUMBATCHES: ", n

        j = 0

        #what does this do?
        L_fn = self.make_L_fn(loss, params)

        print("Starting training...")
        for epoch in range(n_epochs):

            if EXTRA_INFO:
                flog.write("Epoch:{:.2f}\n".format(epoch))

            t = time.time()

            train_err = 0
            train_acc = 0
            train_batches = 0

            # get the batches instances from the train set
            if self.non_uniform_prob:
                batches = self.iterate_minibatches(X_train, Y_train, batch_size)
            else:
                batches = iterate_minibatches(X_train, Y_train, batch_size, shuffle=True)

            m_list = list(range(1,self.m))
            m_tilde = random.choice(m_list)

            for batch in batches:

                # recompute the \mu based on current values of params every m batches
                if j % m_tilde == 0:
                    m_tilde = random.choice(m_list)
                    for mu in self.mu:
                        mu.set_value(0 * mu.get_value())

                    mu_batch = random_minibatches(X_train, Y_train, batch_size, shuffle=False):
                    inputs, targets = mu_batch
                    train_mu(inputs, targets)

                    # ??? where is the summation for mu?
                    # ??? where is the definition of mu? is it a shared theano variabe defined somewhere else?
                    # ??? There is a clean mu to [] in the update_mu, how can it succeed values from other batches?
                    # for mu in self.mu:
                    #     mu.set_value(mu.get_value() / n)
                        # n is number of batches

                j += 1
                # construct the current train sample based on one batch
                inputs, targets = batch
                # print "learning_rate: ", 1. / self.L.get_value()

                # what is L?
                # A:L is learning rate. Ls is the list storing learning rate.
                L = self.Ls[self.idx]
                self.L.set_value(L)

                current_loss, current_acc = val_fn(inputs, targets)

                # Use adaptive option to find the right learning rate by line search
                l_iter = 0
                if self.adaptive:
                    while True:
                        loss_next, sq_sum = L_fn(inputs, targets)
                        if loss_next <= current_loss - 0.5 * sq_sum / self.L.get_value():
                            break
                        else:
                            self.L.set_value(self.L.get_value() * 2)

                        l_iter += 1

                # Print the extra debug infomation
                if EXTRA_INFO:
                    print >>flog, "No. of batch:",train_batches
                    #Each iteration here decrease learning rate by half
                    flog.write("Iteration of L:{:.2f}\n".format(l_iter))
                    print >>flog, "learning_rate: ", 1. / self.L.get_value()

                # Batch updates for parameters w
                train_err += train_w(inputs, targets)
                train_acc += current_acc

                self.Ls[self.idx] = self.L.get_value()
                train_batches += 1

            val_err = 0
            val_acc = 0
            val_batches = 0
            for i, batch in enumerate(iterate_minibatches(X_val, Y_val, batch_size, shuffle=False)):
                inputs, targets = batch
                current_err, current_acc = val_fn(inputs, targets)
                val_err += current_err
                val_acc += current_acc
               # val_err += val_fn(np.array(inputs.todense(), dtype=np.float32), np.array(targets, dtype=np.int32))
                val_batches += 1

            test_err = 0
            test_acc = 0
            test_batches = 0
            for i, batch in enumerate(iterate_minibatches(X_test, y_test, batch_size, shuffle=False)):
                inputs, targets = batch
                current_err, current_acc = val_fn(inputs, targets)
                test_err += current_err
                test_acc += current_acc
#                test_err += test_fn(np.array(inputs.todense(), dtype=np.float32), np.array(targets, dtype=np.int32))
                test_batches += 1

            times.append(time.time() - t)
            print("Epoch {} of {} took {:.3f}s".format(epoch + 1, n_epochs, time.time() - t))

            print("  training loss:\t\t{:.6f}".format(train_err / train_batches))
            print("  validation loss:\t\t{:.6f}".format(val_err / val_batches))
            print("  test loss:\t\t\t{:.6f}".format(test_err / test_batches))
            print("  train accuracy:\t\t{:.6f}".format(train_acc / train_batches))
            print("  validation accuracy:\t\t{:.6f}".format(val_acc / val_batches))
            print("  test accuracy:\t\t{:.6f}\n".format(test_acc / test_batches))

            train_error.append(train_err / train_batches)
            # what is the second parameter here?
            validation_error.append((val_err / val_batches, self.counted_gradient.get_value()))

            acc_train.append(train_acc / train_batches)
            acc_val.append(val_acc / val_batches)
            test_error.append(test_err / test_batches)
            acc_test.append(test_acc / test_batches)

#            if X_val is not None:
#                print("  validation loss:\t\t{:.6f}".format(val_err / val_batches))

        print("Average time per epoch \t {:.3f}".format(np.mean(times)))
        return train_error, validation_error, acc_train, acc_val, acc_test, test_error

    def make_updates(self, loss, params):

        mu_updates = self.make_mu_updates(loss, params)
        w_updates = self.make_w_updates(loss, params)

        return w_updates, mu_updates

    def make_mu_updates(self, loss, params):
        mu_updates = OrderedDict()
        # mu is the partial derivative for every param

        grads = theano.grad(loss, params)
        # grads: the partial gradient w.r.t each param

        self.mu = []
        for param, grad in zip(params, grads):
            # param: for each param in all trainable parameters

            value = param.get_value(borrow=True)

            mu_updates[self.counted_gradient] = self.counted_gradient + 1

            mu = theano.shared(np.zeros(value.shape, dtype=value.dtype), broadcastable=param.broadcastable)
            # seems that mu is just a zero
            # so when we create this mu here???

            mu_updates[mu] = mu + grad
            self.mu.append(mu)

        return mu_updates

    def make_w_updates(self, loss, params):
        w_updates = OrderedDict()

        params_tilde = [theano.shared(x.get_value()) for x in params]
        loss_tilde = theano.clone(loss, replace=zip(params, params_tilde))

        grads = theano.grad(loss, params)
        # grades is the gradient that change per iteration

        grads_tilde = theano.grad(loss_tilde, params_tilde)
        # grades_tilde is the gradient that change every 'm' iterations

        it_num = theano.shared(np.cast['int16'](0))
        # do we succeed the it_num among different update processes??? Otherwise how do we keep counting it correctly?
        it = it_num + 1

        for param, grad, mu, param_tilde, grad_tilde in zip(params, grads, self.mu, params_tilde, grads_tilde):

            new_param = param - (1. / self.L) * (grad - grad_tilde + mu)
            # compute the w_t based on w_{t-1}, w_tilde and mu
            w_updates[param] = new_param
            # update w_t
            w_updates[param_tilde] = ifelse(T.eq(it % self.m, 0), new_param, param_tilde)
            # update w_tilde every 'm' iterations
            # We can print the output of ifelse here to see if we are counting the iterations correctly

            w_updates[self.counted_gradient] = self.counted_gradient + 2
            # why do we have this ???

        if self.adaptive:
            w_updates[self.L] = self.L / 2
            # so why adaptive means we are using larger step sizes?

        self.it_num = it_num

        w_updates[it_num] = it
        return w_updates

    def make_L_fn(self, loss, params):
        grads = theano.grad(loss, params)

        params_next = [x - 1. / self.L * g for x, g in zip(params, grads)]
        loss_next = theano.clone(loss, replace=zip(params, params_next))
        sq_sum = sum((g**2).sum() for g in grads)

        return theano.function([self.input_var, self.target_var], [loss_next, sq_sum])

    def iterate_minibatches(self, inputs, targets, batchsize):
        assert inputs.shape[0] == len(targets)

        indices = np.arange(inputs.shape[0] / batchsize)

        lipschitz_prob = self.get_prob()

        for start_idx in range(0, inputs.shape[0] - batchsize + 1, batchsize):
            [idx] = np.random.choice(indices, size=1, p=lipschitz_prob)
            excerpt = np.arange(inputs.shape[0])[idx * batchsize : (idx + 1) * batchsize]
            self.idx = idx
            yield inputs[excerpt], targets[excerpt]

    def random_minibatches(self, inputs, targets, batchsize):
        assert inputs.shape[0] == len(targets)
        indices = np.arange(inputs.shape[0] / batchsize)
        lipschitz_prob = self.get_prob()
        idx_list = range(0, inputs.shape[0] - batchsize + 1, batchsize)
        while 1:
            start_idx = random.choice(idx_list)
            [idx] = np.random.choice(indices, size=1, p=lipschitz_prob)
            excerpt = np.arange(inputs.shape[0])[idx * batchsize : (idx + 1) * batchsize]
            self.idx = idx
            yield inputs[excerpt], targets[excerpt]

    def get_prob(self):
#        s = sum(L.get_value() for L in self.Ls)
#        return [L.get_value() / s for L in self.Ls]
        s = sum(self.Ls)
        return [L / s for L in self.Ls]

max_stars_repo_name: lizhangzhan/svrg

max_issues_repo_name: lizhangzhan/svrg

max_forks_repo_name: lizhangzhan/svrg

reason: import numpy,import theano,from theano


text:
# This script can go through an instance catalog and apply some post-process
# corrections in order to avoid having to regenerate the whole catalog
import os
import numpy as np
import argparse
import contextlib
from multiprocessing import Pool
import gzip

@contextlib.contextmanager
def fopen(filename, **kwds):
    """
    Return a file descriptor-like object that closes the underlying
    file descriptor when used with the with-statement.
    Parameters
    ----------
    filename: str
        Filename of the instance catalog.
    **kwds: dict
        Keyword arguments to pass to the gzip.open or open functions.
    Returns
    -------
    generator: file descriptor-like generator object that can be iterated
        over to return the lines in a file.
    """
    abspath = os.path.split(os.path.abspath(filename))[0]
    try:
        if filename.endswith('.gz') and os.path.exists(filename):
                fd = gzip.open(filename, **kwds)
        else:
            fd = open(filename.split('.gz')[0], **kwds)
        yield fopen_generator(fd, abspath, **kwds)
    finally:
        fd.close()

def fopen_generator(fd, abspath, **kwds):
    """
    Return a generator for the provided file descriptor that knows how
    to recursively read in instance catalogs specified by the
    includeobj directive.
    """
    with fd as input_:
        for line in input_:
            yield line

def metadata_from_file(file_name):
    """
    Read in the InstanceCatalog specified by file_name.
    Return a dict of the header values from that InstanceCatalog.
    Simpler version than the ImSim code
    """
    input_params = {}
    catalog_files = []
    with fopen(file_name, mode='rt') as in_file:
        for line in in_file:
            if line[0] == '#':
                continue

            params = line.strip().split()

            if params[0] == 'includeobj':
                catalog_files.append(params[1])
                continue

            float_val = float(params[1])
            int_val = int(float_val)
            if np.abs(float_val-int_val) > 1.0e-10:
                val = float_val
            else:
                val = int_val
            input_params[params[0]] = val

    commands = dict(((key, value) for key, value in input_params.items()))
    return commands, catalog_files

def apply_extinction_correction(tokens):
    """
    This function updates the extinction in the provided entry
    """
    # Extracting the extinction parameters:
    if tokens[17].lower() != 'none':
        internal_av = float(tokens[18])
        internal_rv = float(tokens[19])
    else:
        internal_av = 0
        internal_rv = 0

    # Proposed staged cut
    corrected = False
    if internal_rv < 0.1:
        internal_rv = 0.1
        internal_av = np.clip(internal_av, 0.0, 1.0)
        corrected = True
    elif internal_rv < 1:
        if internal_av > 1:
            corrected = True
        internal_av = np.clip(internal_av, 0.0, 1.0)
    elif internal_av < 0:
        # I'm not counting the lower cut on Av, almost without effect for bulges
        # and likely to be minimal for disks
        internal_av = 0

    # update tokens
    if corrected:
        tokens[18] = ("%.9f" % internal_av).rstrip('0')
        tokens[19] = ("%.9f" % internal_rv).rstrip('0')

    return corrected, tokens

def fix_disk_knots(in_instcat_disk, in_instcat_knots,
                   out_instcat_disk, out_instcat_knots):

    # Use .fopen to read in the command and object lines from the
    # instance catalog.
    count_extinction = 0
    count_line = 0
    with fopen(in_instcat_disk, mode='rt') as input_disk,   \
         fopen(in_instcat_knots, mode='rt') as input_knots, \
         open(out_instcat_disk, 'w') as output_disk, \
         open(out_instcat_knots, 'w') as output_knots:

        # We first go through the knots catalog because some entries are missing
        # compared to the full bulge/disk catalog (faint knots have already been
        # removed at the instance catalog creation level)
        for line_knots in input_knots:

            # Extract the galaxy ID for that knots component
            tokens_knots = line_knots.strip().split()
            id_knots = int(tokens_knots[1]) >> 10

            found = False
            # Loop through the disk catalogs
            for line_disk in input_disk:
                tokens_disk = line_disk.strip().split()
                id_disk = int(tokens_disk[1]) >> 10

                # If the galaxy is offensive, clip the av and rv values
                corrected, tokens_disk = apply_extinction_correction(tokens_disk)
                # If the galaxy is offensive, clip the av and rv values
                if corrected:
                    count_extinction += 1

                count_line +=1
                if id_disk == id_knots:
                    found=True
                    break
                else:
                    line_disk = ' '.join(tokens_disk)
                    output_disk.write(line_disk.strip()+'\n')

            if not found:
                print("ERROR: object ids do not match between input knots and disks catalogs")
                exit(-1)

            # Get total flux
            magnorm_disk = np.float(tokens_disk[4])
            magnorm_knots = np.float(tokens_knots[4])
            total_flux = 10.**(-magnorm_disk/2.5) + 10.**(-magnorm_knots/2.5)
            knots_flux_ratio = 10.**(-magnorm_knots/2.5) / total_flux

            # Apply smooth flux ration scaling
            size = np.float(tokens_disk[13])
            knots_flux_ratio = knots_flux_ratio * 0.5*(1 - np.tanh(np.log(size/1.5)))

            magnorm_disk = -2.5*np.log10((1-knots_flux_ratio)*total_flux)
            magnorm_knots = -2.5*np.log10(knots_flux_ratio*total_flux)

            # Update the entry
            tokens_disk[4] = ("%.9f" % magnorm_disk).rstrip('0')
            tokens_knots[4] = ("%.9f" % magnorm_knots).rstrip('0')
            # Making sure that the extinction paramters remain the same between disk and knots
            tokens_knots[18] = tokens_disk[18]
            tokens_knots[19] = tokens_disk[19]
            line_disk = ' '.join(tokens_disk)
            line_knots = ' '.join(tokens_knots)

            # Write the catalogs
            output_disk.write(line_disk.strip()+'\n')
            output_knots.write(line_knots.strip()+'\n')

    print("Fixed extinction for %d disks bulge out of %d"%(count_extinction, count_line))

def fix_bulge(in_instcat_bulge, out_instcat_bulge):
    count_extinction = 0
    count_line =0
    with fopen(in_instcat_bulge, mode='rt') as input_bulge,  \
         open(out_instcat_bulge, 'w') as output_bulge:

        # Now we go through the bulge catalog independently, fixing the extinction
        for line_bulge in input_bulge:
            tokens_bulge = line_bulge.strip().split()
            id_bulge = int(tokens_bulge[1]) >> 10

            corrected, tokens_bulge = apply_extinction_correction(tokens_bulge)
            if corrected:
                count_extinction += 1
            count_line +=1
            line_bulge = ' '.join(tokens_bulge)
            output_bulge.write(line_bulge.strip()+'\n')
    print("Fixed extinction for %d bulge out of %d"%(count_extinction, count_line))

def process_instance_catalog(args):
    """
    Function that processes a single instance catalog
    """
    input_cat, output_path = args
    # Find the visit id
    metadata, catalog_files = metadata_from_file(input_cat)
    visitID = metadata['obshistid']
    input_path = input_cat.split('/')[:-1]
    input_path = '/'.join(input_path)
    print('Copying catalog from %s'%input_path)

    # Create output directory
    output_path = os.path.join(output_path,'{0:08d}'.format(int(visitID)))
    if not os.path.exists(output_path):
        os.makedirs(output_path)

    # Copy over the content of the instance catalog
    os.system("cp -ar %s/* %s"%(input_path, output_path))

    # Checking that all the files are gzipped
    for f in catalog_files:
        basename = f.split('.gz')[0]
        # Skipping these files as we will gzip them at the end anyway
        if 'disk' in basename or 'bulge' in basename or 'knots' in basename:
            continue
        name = os.path.join(output_path,basename)
        if os.path.exists(name):
            print('Gzipping '+ name)
            os.system("gzip -f %s " % name)

    # Processes catalogs, prefer the ungzipped version if it exists
    input_disk=output_path+'/disk_gal_cat_%d.txt'%visitID
    input_bulge=output_path+'/bulge_gal_cat_%d.txt'%visitID
    input_knots=output_path+'/knots_cat_%d.txt'%visitID

    # Checking that the gz files exist, otherwise remove the gz extension
    if not os.path.exists(input_disk):
        input_disk=output_path+'/disk_gal_cat_%d.txt.gz'%visitID
    if not os.path.exists(input_bulge):
        input_bulge=output_path+'/bulge_gal_cat_%d.txt.gz'%visitID
    if not os.path.exists(input_knots):
        input_knots=output_path+'/knots_cat_%d.txt.gz'%visitID

    output_disk=output_path+'/disk_gal_cat_%d.txt'%visitID
    output_bulge=output_path+'/bulge_gal_cat_%d.txt'%visitID
    output_knots=output_path+'/knots_cat_%d.txt'%visitID

    tmp_disk=output_path+'/tmp_disk_gal_cat_%d.txt'%visitID
    tmp_bulge=output_path+'/tmp_bulge_gal_cat_%d.txt'%visitID
    tmp_knots=output_path+'/tmp_knots_cat_%d.txt'%visitID

    print('Processing disks and knots for %d'%visitID)
    fix_disk_knots(input_disk, input_knots, tmp_disk, tmp_knots)
    os.system("mv %s %s"%(tmp_disk, output_disk))
    os.system("mv %s %s"%(tmp_knots, output_knots))

    print('Processing bulges')
    fix_bulge(input_bulge, tmp_bulge)
    os.system("mv %s %s"%(tmp_bulge, output_bulge))
    print('Gzipping....')
    os.system("gzip -f %s %s %s" % (output_disk, output_bulge, output_knots))

    print('Done.')

if __name__ == '__main__':
    parser = argparse.ArgumentParser(description='Instance catalog crawler applying corrections in post-processing')
    parser.add_argument('input_cats', type=str,help='List of instance catalogs')
    parser.add_argument('output_path', type=str, help='Directory in which to store the corrected catalog')
    args = parser.parse_args()

    filenames = []
    with open(args.input_cats, 'r') as f:
        for line in f:
            if len(line.strip()) > 0:
                filenames.append(line.strip())

    p = Pool(24)
    p.map(process_instance_catalog, [(f,args.output_path) for f in filenames] )


text:
import numpy as np
from flare.modules import otf_parser, analyze_md
from flare import kernels, otf, env
import concurrent.futures
import time

class TestOTF:
    def __init__(self, otf_file, gp_cell, kernel,
                 kernel_grad, cutoffs, start_list, test_snaps):
        self.test_snaps = test_snaps
        self.struc = None  # populate later

        # make gp and aimd run
        otf_run = otf_parser.OtfAnalysis(otf_file)
        call_no = len(otf_run.gp_position_list)
        self.gp_model = otf_run.make_gp(gp_cell, kernel, kernel_grad, 'BFGS',
                                        call_no, start_list, cutoffs)

    # adapted from analyze_gp.py
    def predict_forces_on_test_set(self, aimd_file, aimd_cell):
        # get aimd trajectory
        aimd_run = analyze_md.MDAnalysis(aimd_file, aimd_cell)

        all_predictions = np.array([])
        all_stds = np.array([])
        all_forces = np.array([])

        for snap in self.test_snaps:
            # get structure and forces from AIMD trajectory
            self.struc = aimd_run.get_structure_from_snap(snap)
            forces_curr = aimd_run.get_forces_from_snap(snap)

            # predict forces and stds
            self.predict_on_structure_par()
            print(time.time())
            predictions = self.struc.forces
            stds = self.struc.stds

            # append results
            all_predictions = np.append(all_predictions, predictions)
            all_stds = np.append(all_stds, stds)
            all_forces = np.append(all_forces, forces_curr)

        return all_predictions, all_stds, all_forces

    def predict_on_structure_par(self):
        atom_list = list(range(self.struc.positions.shape[0]))
        n = 0
        with concurrent.futures.ProcessPoolExecutor() as executor:
            for res in executor.map(self.predict_on_atom, atom_list):
                for i in range(3):
                    self.struc.forces[n][i] = res[0][i]
                    self.struc.stds[n][i] = res[1][i]
                n += 1

    def predict_on_structure(self):
        for n in range(self.struc.nat):
            chemenv = \
                env.AtomicEnvironment(self.struc, n, self.gp_model.cutoffs)
            for i in range(3):
                force, var = self.gp_model.predict(chemenv, i + 1)
                self.struc.forces[n][i] = float(force)
                self.struc.stds[n][i] = np.sqrt(np.absolute(var))

    # take prediction functions from otf
    def predict_on_atom(self, atom):
        chemenv = \
            env.AtomicEnvironment(self.struc, atom, self.gp_model.cutoffs)
        comps = []
        stds = []
        # predict force components and standard deviations
        for i in range(3):
            force, var = self.gp_model.predict(chemenv, i+1)
            comps.append(float(force))
            stds.append(np.sqrt(np.absolute(var)))

        return comps, stds


text:
import numpy as np
from PIL import Image

import torch
class LFW(object):
    def __init__(self, imgl, imgr):

        self.imgl_list = imgl
        self.imgr_list = imgr

    def __getitem__(self, index):
        imgl = np.array(Image.open(self.imgl_list[index]))
        if len(imgl.shape) == 2:
            imgl = np.stack([imgl] * 3, 2)
        imgr = np.array(Image.open(self.imgr_list[index]))
        if len(imgr.shape) == 2:
            imgr = np.stack([imgr] * 3, 2)

        # imgl = imgl[:, :, ::-1]
        # imgr = imgr[:, :, ::-1]
        imglist = [imgl, imgl[:, ::-1, :], imgr, imgr[:, ::-1, :]]
        for i in range(len(imglist)):
            imglist[i] = (imglist[i] - 127.5) / 128.0
            imglist[i] = imglist[i].transpose(2, 0, 1)
        imgs = [torch.from_numpy(i).float() for i in imglist]
        return imgs

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


text:
#Some notebook setup stuff

#we may need some code in the ../python directory and/or matplotlib styles
import sys
import os
sys.path.append('../python/')

#set up matplotlib
os.environ['MPLCONFIGDIR'] = '../mplstyles'
#print(os.environ['MPLCONFIGDIR'])
import matplotlib as mpl
from matplotlib import pyplot as plt
#got smarter about the mpl config: see mplstyles/ directory
plt.style.use('standard')
#print(mpl.__version__)
#print(mpl.get_configdir())

#fonts
# Set the font dictionaries (for plot title and axis titles)
title_font = {'fontname':'Arial', 'size':'16', 'color':'black', 'weight':'normal',
              'verticalalignment':'bottom'} # Bottom vertical alignment for more space
axis_font = {'fontname':'Arial', 'size':'32'}
legend_font = {'fontname':'Arial', 'size':'22'}

#fonts global settings
mpl.rc('font',family=legend_font['fontname'])

#set up numpy
import numpy as np
import pandas as pd

max_stars_repo_name: villano-lab/nrSiCap

max_issues_repo_name: villano-lab/nrSiCap

max_forks_repo_name: villano-lab/nrSiCap