text:
# -*- coding: utf-8 -*-
# file: model.py
# author: JinTian
# time: 07/03/2017 3:07 PM
# Copyright 2017 JinTian. All Rights Reserved.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
# ------------------------------------------------------------------------
import tensorflow as tf
import numpy as np

def rnn_model(model, input_data, output_data, vocab_size, rnn_size=128, num_layers=2, batch_size=64,
              input_num=128,output_num=1,learning_rate=0.01,use_cnn=True):
    """
    construct rnn seq2seq model.
    :param model: model class
    :param input_data: input data placeholder
    :param output_data: output data placeholder
    :param vocab_size:
    :param rnn_size:
    :param num_layers:
    :param batch_size:
    :param learning_rate:
    :return:
    """
    end_points = {}

    if model == 'rnn':
        cell_fun = tf.contrib.rnn.BasicRNNCell
        cell = cell_fun(rnn_size)
    elif model == 'gru':
        cell_fun = tf.contrib.rnn.GRUCell
        cell = cell_fun(rnn_size)
    elif model == 'lstm':
        cell_fun = tf.contrib.rnn.BasicLSTMCell
        cell = cell_fun(rnn_size, state_is_tuple=True)

    cell = tf.contrib.rnn.MultiRNNCell([cell] * num_layers, state_is_tuple=True)

    if output_data is not None:
        initial_state = cell.zero_state(batch_size, tf.float32)
    else:
        initial_state = cell.zero_state(1, tf.float32)
    if use_cnn:
        with tf.name_scope('fc1'):
            fc1_weights = tf.Variable(  # fully connected, depth 512.
                tf.truncated_normal([input_num, 128],
                                    # mean=1.0,
                                    stddev=0.1,
                                    dtype=tf.float32))
            fc1_biases = tf.Variable(tf.constant(1.0, shape=[128], dtype=tf.float32))
            fc1 = tf.nn.relu(tf.matmul(tf.to_float(input_data), fc1_weights) + fc1_biases)
        with tf.name_scope('fc2'):
            fc2_weights = tf.Variable(  # fully connected, depth 512.
                tf.truncated_normal([128, 256],
                                    # mean=1.0,
                                    stddev=0.1,
                                    dtype=tf.float32))
            fc2_biases = tf.Variable(tf.constant(1.0, shape=[256], dtype=tf.float32))
            fc2 = tf.nn.relu(tf.matmul(tf.to_float(fc1), fc2_weights) + fc2_biases)
        with tf.name_scope('fc3'):
            fc3_weights = tf.Variable(  # fully connected, depth 512.
                tf.truncated_normal([256, 512],
                                    # mean=1.0,
                                    stddev=0.1,
                                    dtype=tf.float32))
            fc3_biases = tf.Variable(tf.constant(1.0, shape=[512], dtype=tf.float32))
            fc3 = tf.nn.relu(tf.matmul(tf.to_float(fc2), fc3_weights) + fc3_biases)
        with tf.name_scope('fc4'):
            fc4_weights = tf.Variable(  # fully connected, depth 512.
                tf.truncated_normal([512, 128*output_num],
                                    # mean=1.0,
                                    stddev=0.1,
                                    dtype=tf.float32))
            fc4_biases = tf.Variable(tf.constant(1.0, shape=[128*output_num], dtype=tf.float32))
            fc4 = tf.nn.relu(tf.matmul(fc3, fc4_weights) + fc4_biases)

        # [batch_size, ?, rnn_size] = [64, ?, 128]
        embedded=tf.reshape(fc4,[batch_size,-1,128])
        outputs, last_state = tf.nn.dynamic_rnn(cell, embedded, initial_state=initial_state)
    else:
        with tf.device("/cpu:0"):
            embedding = tf.get_variable('embedding', initializer=tf.random_uniform(
                [vocab_size + 1, rnn_size], -1.0, 1.0))
            inputs = tf.nn.embedding_lookup(embedding, input_data)
    # embedded = tf.reshape(input_data, [batch_size, -1, 128])
    # outputs, last_state = tf.nn.dynamic_rnn(cell, embedded, initial_state=initial_state)
    # outputs, last_state = tf.nn.dynamic_rnn(cell, inputs, initial_state=initial_state)
    output = tf.reshape(outputs, [-1, rnn_size])

    weights = tf.Variable(tf.truncated_normal([rnn_size, vocab_size + 1]))
    bias = tf.Variable(tf.zeros(shape=[vocab_size + 1]))
    logits = tf.nn.bias_add(tf.matmul(output, weights), bias=bias)
    # [?, vocab_size+1]

    if output_data is not None:
        # output_data must be one-hot encode
        labels = tf.one_hot(tf.reshape(output_data, [-1]), depth=vocab_size + 1)
        # should be [?, vocab_size+1]
        loss = tf.nn.softmax_cross_entropy_with_logits(labels=labels, logits=logits)

        # loss shape should be [?, vocab_size+1]
        total_loss = tf.reduce_mean(loss)
        regularizers = (tf.nn.l2_loss(fc1_weights) + tf.nn.l2_loss(fc1_biases) +
                tf.nn.l2_loss(fc2_weights) + tf.nn.l2_loss(fc2_biases)+
                tf.nn.l2_loss(fc3_weights) + tf.nn.l2_loss(fc3_biases)+
                tf.nn.l2_loss(fc4_weights) + tf.nn.l2_loss(fc4_biases))
        total_loss+=5e-4 * regularizers
        train_op = tf.train.AdamOptimizer(learning_rate).minimize(total_loss)

        end_points['initial_state'] = initial_state
        end_points['output'] = output
        end_points['train_op'] = train_op
        end_points['total_loss'] = total_loss
        end_points['loss'] = loss
        end_points['last_state'] = last_state
        prediction = tf.nn.softmax(logits)
        end_points['prediction'] = prediction
    else:
        prediction = tf.nn.softmax(logits)

        end_points['initial_state'] = initial_state
        end_points['last_state'] = last_state
        end_points['prediction'] = prediction

    return end_points

max_stars_repo_name: gnosis-xian/ssq


text:
# -*- coding: utf-8 -*-

# Copyright (C) 2017-2018 Andreas Franck and Giacomo Costantini
# Copyright (C) 2017-2018 University of Southampton

# VISR Binaural Synthesis Toolkit (BST)
# Authors: Andreas Franck and Giacomo Costantini
# Project page: http://cvssp.org/data/s3a/public/BinauralSynthesisToolkit/

# The Binaural Synthesis Toolkit is provided under the ISC (Internet Systems Consortium) license
# https://www.isc.org/downloads/software-support-policy/isc-license/ :

# Permission to use, copy, modify, and/or distribute this software for any
# purpose with or without fee is hereby granted, provided that the above
# copyright notice and this permission notice appear in all copies.
#
# THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
# WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY
# AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT,
# INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS
# OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS
# ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.

# We kindly ask to acknowledge the use of this software in publications or software.
# Paper citation:
# Andreas Franck, Giacomo Costantini, Chris Pike, and Filippo Maria Fazi,
# “An Open Realtime Binaural Synthesis Toolkit for Audio Research,” in Proc. Audio Eng.
# Soc. 144th Conv., Milano, Italy, 2018, Engineering Brief.
# http://www.aes.org/e-lib/browse.cfm?elib=19525

# The Binaural Synthesis Toolkit is based on the VISR framework. Information about the VISR,
# including download, setup and usage instructions, can be found on the VISR project page
# http://cvssp.org/data/s3a/public/VISR .

import numpy as np

# Core VISR packages
import visr
import rbbl
import pml
import rcl

from visr_bst.hoa_components import HoaObjectEncoder
from visr_bst.hoa_components import HoaCoefficientRotation
from visr_bst.util import readSofaFile

class HoaObjectToBinauralRenderer( visr.CompositeComponent ):
    """
    Component to render binaural audio from plane wave and point source objects using an Higher Order Ambisonics (HOA)
    algorithm.
    """
    def __init__( self,
                 context, name, parent,
                 numberOfObjects,
                 maxHoaOrder = None,
                 sofaFile = None,
                 decodingFilters = None,
                 interpolationSteps = None,
                 headOrientation = None,
                 headTracking = True,
                 objectChannelAllocation = False,
                 fftImplementation = 'default'
                 ):
        """
        Constructor.

        Parameters
        ----------
        context : visr.SignalFlowContext
            Standard visr.Component construction argument, holds the block size and the sampling frequency
        name : string
            Name of the component, Standard visr.Component construction argument
        parent : visr.CompositeComponent
            Containing component if there is one, None if this is a top-level component of the signal flow.
        numberOfObjects : int
            The number of audio objects to be rendered.
        maxHoaOrder: int or None
            The maximum HOA order that can be reproduced. If None, the HOA order is deduced
            from the first dimension of the HOA filters (possibly contained in a SOFA file).
        sofaFile: string or NoneType
        decodingFilters : numpy.ndarray or NoneType
            Alternative way to provide the HOA decoding filters.
        interpolationSteps: int
        headOrientation : array-like
            Head orientation in spherical coordinates (2- or 3-element vector or list). Either a static orientation (when no tracking is used),
            or the initial view direction
        headTracking: bool
            Whether dynamic head tracking is active.
        objectChannelAllocation: bool
            Whether the processing resources are allocated from a pool of resources
            (True), or whether fixed processing resources statically tied to the audio signal channels are used.
            Not implemented at the moment, so leave the default value (False).
        fftImplementation: string, optional
            The FFT library to be used in the filtering. THe default uses VISR's
            default implementation for the present platform.
        """
        if (decodingFilters is None) == (sofaFile is None ):
            raise ValueError( "HoaObjectToBinauralRenderer: Either 'decodingFilters' or 'sofaFile' must be provided." )
        if sofaFile is None:
            filters = decodingFilters
        else:
            # pos and delays are not used here.
            [pos, filters, delays] = readSofaFile( sofaFile )

        if maxHoaOrder is None:
            numHoaCoeffs = filters.shape[0]
            orderP1 = int(np.floor(np.sqrt(numHoaCoeffs)))
            if orderP1**2 != numHoaCoeffs:
                raise ValueError( "If maxHoaOrder is not given, the number of HOA filters must be a square number" )
            maxHoaOrder = orderP1 - 1
        else:
            numHoaCoeffs = (maxHoaOrder+1)**2

        if filters.ndim != 3 or filters.shape[1] != 2 or filters.shape[0] < numHoaCoeffs:
            raise ValueError( "HoaObjectToBinauralRenderer: the filter data must be a 3D matrix where the second dimension is 2 and the first dimension is equal or larger than (maxHoaOrder+1)^2." )

        super( HoaObjectToBinauralRenderer, self ).__init__( context, name, parent )
        self.objectSignalInput = visr.AudioInputFloat( "audioIn", self, numberOfObjects )
        self.binauralOutput = visr.AudioOutputFloat( "audioOut", self, 2 )
        self.objectVectorInput = visr.ParameterInput( "objectVector", self, pml.ObjectVector.staticType,
                                                     pml.DoubleBufferingProtocol.staticType,
                                                     pml.EmptyParameterConfig() )

        if interpolationSteps is None:
            interpolationSteps = context.period

        self.objectEncoder = HoaObjectEncoder( context, 'HoaEncoder', self,
                                              numberOfObjects = numberOfObjects,
                                              hoaOrder = maxHoaOrder,
                                              channelAllocation = objectChannelAllocation )

        self.parameterConnection( self.objectVectorInput, self.objectEncoder.parameterPort("objectVector"))

        self.encoderMatrix = rcl.GainMatrix( context, "encoderMatrix", self,
                                             numberOfInputs = numberOfObjects,
                                             numberOfOutputs = (maxHoaOrder+1)**2,
                                             interpolationSteps = interpolationSteps,
                                             initialGains = 0.0,
                                             controlInput = True )
        self.audioConnection( self.objectSignalInput, self.encoderMatrix.audioPort("in") )

        filterMtx = np.concatenate( (filters[0:numHoaCoeffs,0,:], filters[0:numHoaCoeffs,1,:]) )

        routings = rbbl.FilterRoutingList()
        for idx in range(0,numHoaCoeffs):
            routings.addRouting( idx, 0, idx, 1.0 )
            routings.addRouting( idx, 1, idx+numHoaCoeffs, 1.0 )

        self.binauralFilterBank = rcl.FirFilterMatrix( context, 'binauralFilterBank', self,
                                                    numberOfInputs = numHoaCoeffs,
                                                    numberOfOutputs = 2,
                                                    filterLength = filters.shape[-1],
                                                    maxFilters = 2*numHoaCoeffs,
                                                    maxRoutings = 2*numHoaCoeffs,
                                                    filters = filterMtx,
                                                    routings = routings,
                                                    controlInputs=rcl.FirFilterMatrix.ControlPortConfig.NoInputs,
                                                    fftImplementation = fftImplementation )

        self.audioConnection( self.encoderMatrix.audioPort("out"), self.binauralFilterBank.audioPort("in") )
        self.audioConnection( self.binauralFilterBank.audioPort("out"), self.binauralOutput )

        if headTracking:
            self.trackingInput = visr.ParameterInput( "tracking", self, pml.ListenerPosition.staticType,
                                              pml.DoubleBufferingProtocol.staticType,
                                              pml.EmptyParameterConfig() )
            self.coefficientRotator = HoaCoefficientRotation( context, 'coefficientRotator', self,
                                                              numberOfObjects = numberOfObjects,
                                                              hoaOrder = maxHoaOrder)
            self.parameterConnection( self.trackingInput, self.coefficientRotator.parameterPort("tracking") )
            self.parameterConnection( self.objectEncoder.parameterPort("coefficientOutput"),
                                     self.coefficientRotator.parameterPort("coefficientInput") )
            self.parameterConnection( self.coefficientRotator.parameterPort("coefficientOutput"),
                                     self.encoderMatrix.parameterPort("gainInput") )
        else:
            self.parameterConnection( self.objectEncoder.parameterPort("coefficientOutput"),
                                     self.encoderMatrix.parameterPort("gainInput") )

ext: py

reason: import numpy


text:
#!/usr/bin/env python3.5
"""Test ModelPart class."""

import os
import tempfile
import unittest

import numpy as np
from numpy.testing import assert_array_equal
import tensorflow as tf

from neuralmonkey.vocabulary import Vocabulary
from neuralmonkey.encoders.recurrent import SentenceEncoder
from neuralmonkey.model.sequence import EmbeddedSequence

class Test(unittest.TestCase):
    """Test capabilities of model part."""

    @classmethod
    def setUpClass(cls):
        tf.reset_default_graph()
        cls.dataset = {
            "id": tf.constant([["hello", "world"], ["test", "this"]]),
            "data_id": tf.constant([["A", "B", "C"], ["D", "E", "F"]])}

    def test_reuse(self):
        vocabulary = Vocabulary(["a", "b"])

        seq1 = EmbeddedSequence(
            name="seq1",
            vocabulary=vocabulary,
            data_id="id",
            embedding_size=10)
        seq1.register_input(self.dataset)

        seq2 = EmbeddedSequence(
            name="seq2",
            vocabulary=vocabulary,
            embedding_size=10,
            data_id="id")
        seq2.register_input(self.dataset)

        seq3 = EmbeddedSequence(
            name="seq3",
            vocabulary=vocabulary,
            data_id="id",
            embedding_size=10,
            reuse=seq1)
        seq3.register_input(self.dataset)

        # blessing
        self.assertIsNotNone(seq1.embedding_matrix)
        self.assertIsNotNone(seq2.embedding_matrix)
        self.assertIsNotNone(seq3.embedding_matrix)

        sess = tf.Session()
        sess.run(tf.global_variables_initializer())

        params = sess.run((seq1.embedding_matrix, seq2.embedding_matrix,
                           seq3.embedding_matrix))

        with self.assertRaises(AssertionError):
            assert_array_equal(params[0], params[1])

        assert_array_equal(params[0], params[2])

    def test_save_and_load(self):
        """Try to save and load encoder."""
        vocabulary = Vocabulary(["a", "b"])

        checkpoint_file = tempfile.NamedTemporaryFile(delete=False)
        checkpoint_file.close()

        encoder = SentenceEncoder(
            name="enc", vocabulary=vocabulary, data_id="data_id",
            embedding_size=10, rnn_size=20, max_input_len=30,
            save_checkpoint=checkpoint_file.name,
            load_checkpoint=checkpoint_file.name)

        encoder.input_sequence.register_input(self.dataset)

        # NOTE: This assert needs to be here otherwise the model has
        # no parameters since the sentence encoder is initialized lazily
        self.assertIsInstance(encoder.temporal_states, tf.Tensor)

        encoders_variables = tf.get_collection(
            tf.GraphKeys.GLOBAL_VARIABLES, scope="enc")

        sess_1 = tf.Session()
        sess_1.run(tf.global_variables_initializer())
        encoder.save(sess_1)

        sess_2 = tf.Session()
        sess_2.run(tf.global_variables_initializer())
        encoder.load(sess_2)

        values_in_sess_1 = sess_1.run(encoders_variables)
        values_in_sess_2 = sess_2.run(encoders_variables)

        self.assertTrue(
            all(np.all(v1 == v2) for v1, v2 in
                zip(values_in_sess_1, values_in_sess_2)))

        os.remove(checkpoint_file.name)

if __name__ == "__main__":
    unittest.main()

max_stars_repo_name: ufal/neuralmonkey

max_issues_repo_name: ufal/neuralmonkey

max_forks_repo_name: ufal/neuralmonkey

reason: import numpy,from numpy


text:
import numpy as np
import pandas as pd
import cellranger.rna.library as rna_library

HIGH_UMI_CORRECTION_THRESHOLD = 0.5
NUM_READS_THRESHOLD = 10000

def detect_aggregate_bcs(correction_data):
    """ Given a csv containing the number of reads and umi corrected reads for each barcode,
        filter out barcodes that exceed pre-defined thresholds of fraction corrected reads and fraction total reads
    """

    summary = correction_data[correction_data['library_type'] == rna_library.ANTIBODY_LIBRARY_TYPE]
    corrected_ratio = summary['umi_corrected_reads'] / summary['reads']
    total_reads = np.float(np.sum(summary['reads']))
    reads_ratio = summary['reads'] / total_reads
    correction_metrics = pd.DataFrame({'Fraction Corrected Reads': corrected_ratio, 'Fraction Total Reads': reads_ratio})
    augmented = summary.join(correction_metrics)

    ### filter out the aggregate barcodes according to set thresholds
    bcs_to_remove = {}
    for idx, row in augmented.iterrows():
        if row['Fraction Corrected Reads'] > HIGH_UMI_CORRECTION_THRESHOLD and row['reads'] > NUM_READS_THRESHOLD:
            bcs_to_remove[row['barcode']] = row['Fraction Total Reads']

    ### make a slice of the barcode summary csv containing just the aggregate barcodes
    augmented.set_index("barcode", inplace=True)
    removed_bcs_df = augmented.loc[bcs_to_remove.keys()].drop(['umis', 'candidate_dup_reads'], axis=1)
    removed_bcs_df.round = removed_bcs_df.round({'Fraction Corrected Reads': 3, 'Fraction Total Reads': 3})

    return bcs_to_remove.keys(), np.sum(bcs_to_remove.values()), removed_bcs_df

max_stars_repo_name: qiangli/cellranger

max_issues_repo_name: qiangli/cellranger

max_forks_repo_name: qiangli/cellranger


text:
# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.11.3
#   kernelspec:
#     display_name: Python 3
#     name: python3
# ---

# + [markdown] id="view-in-github" colab_type="text"
#

# + [markdown] id="rF208fIxvq8m"
# # Introduction to neural networks using Flax
#
#
#
# Flax / Linen is a neural net library, built on top of JAX, "designed to offer an implicit variable management API to save the user from having to manually thread thousands of variables through a complex tree of functions." To handle both current and future JAX transforms (configured and composed in any way), Linen Modules are defined as explicit functions of the form
# $$
# f(v_{in}, x) \rightarrow v_{out}, y
# $$
# Where $v_{in}$ is the collection of variables (eg. parameters) and PRNG state used by the model, $v_{out}$ the mutated output variable collections, $x$ the input data and $y$ the output data. We illustrate this below. Our tutorial is based on the official [flax intro](https://flax.readthedocs.io/en/latest/notebooks/flax_basics.html) and [linen colab](https://github.com/google/flax/blob/master/docs/notebooks/linen_intro.ipynb). Details are in the [flax source code](https://flax.readthedocs.io/en/latest/_modules/index.html). Note: please be sure to read our [JAX tutorial](https://github.com/probml/pyprobml/blob/master/book1/intro/jax_intro.ipynb) first.
#

# + id="uRAzAXYXvztz"
import numpy as np
#np.set_printoptions(precision=3)
np.set_printoptions(formatter={'float': lambda x: "{0:0.5f}".format(x)})

import matplotlib.pyplot as plt

# + colab={"base_uri": "https://localhost:8080/"} id="1ob8P9ALvkcM" outputId="573415d3-eba0-4d6a-e51a-674aea067217"
# Install the latest JAXlib version.
# #!pip install --upgrade -q pip jax jaxlib

# + id="68kI74E1vvEI" colab={"base_uri": "https://localhost:8080/"} outputId="4525fc58-0c43-4909-dad9-1172ef516cbf"
import jax
from jax import lax, random, numpy as jnp
key = random.PRNGKey(0)

# + id="N3dXu6XY6U0H"
from typing import Any, Callable, Dict, Iterator, Mapping, Optional, Sequence, Tuple

# Useful type aliases
Array = jnp.ndarray
PRNGKey = Array
Batch = Mapping[str, np.ndarray]
OptState = Any

# + colab={"base_uri": "https://localhost:8080/"} id="7pcNcE9_Qj_l" outputId="107ed9cb-1d4b-46ab-8032-60836e4b083e"
# Install Flax at head:
# !pip install --upgrade -q git+https://github.com/google/flax.git

# + id="s80k9sonQfDi"
import flax
from flax.core import freeze, unfreeze
from flax import linen as nn
from flax import optim

from jax.config import config
config.enable_omnistaging() # Linen requires enabling omnistaging

# + [markdown] id="aGrUFJYxjyL7"
# # MLP in vanilla JAX
#
# We construct a simple MLP with L hidden layers (relu activation), and scalar output (linear activation).
#
# Note: JAX and Flax, like NumPy, are row-based systems, meaning that vectors are represented as row vectors and not column vectors.
#

# + id="mWQGVJMP0VMB"
# We define the parameter initializers using a signature that is flax-compatible
# https://flax.readthedocs.io/en/latest/_modules/jax/_src/nn/initializers.html

def weights_init(key, shape, dtype=jnp.float32):
  return random.normal(key, shape, dtype)
  #return jnp.ones(shape, dtype)

def bias_init(key, shape, dtype=jnp.float32):
  return jnp.zeros(shape, dtype)

def relu(a):
  return jnp.maximum(a, 0)

# + id="GepkhhTh-9b-"
# A minimal MLP class

class MLP0():
  features: Sequence[int] # number of features in each layer

  def __init__(self, features): # class constructor
    self.features = features

  def init(self, key, x): # initialize parameters
    in_size = np.shape(x)[1]
    sizes = np.concatenate( ([in_size], self.features) )
    nlayers = len(sizes)
    params = {}
    for i in range(nlayers-1):
      in_size = sizes[i]
      out_size = sizes[i+1]
      subkey1, subkey2, key = random.split(key, num=3)
      W = weights_init(subkey1, (in_size, out_size) )
      b = bias_init(subkey2, out_size)
      params[f'W{i}'] = W
      params[f'b{i}'] = b
    return params

  def apply(self, params, x): # forwards pass
    activations = x
    nhidden_layers = len(self.features)-1
    for i in range(nhidden_layers):
      W = params[f'W{i}'];
      b = params[f'b{i}'];
      outputs = jnp.dot(activations, W) + b
      activations = relu(outputs)
    # for final layer, no activation function
    i = nhidden_layers
    outputs = jnp.dot(activations, params[f'W{i}']) + params[f'b{i}']
    return outputs

# + colab={"base_uri": "https://localhost:8080/"} id="8jFS4SNO0V_I" outputId="72d9d449-742e-45d1-9641-a4f3a6390e26"
key = random.PRNGKey(0)
D = 3
N = 2
x = random.normal(key, (N,D,))
layer_sizes = [3,1] # 1 hidden layer of size 3, 1 scalar output

model0 = MLP0(layer_sizes)
params0 = model0.init(key, x)

print('params')
for k,v in params0.items():
  print(k, v.shape)
  print(v)

y0 = model0.apply(params0, x)
print('\noutput')
print(y0)

# + [markdown] id="rBtPT-drBkGA"
# # Our first flax model
#
# Here we recreate the vanilla model in flax. Since we don't specify how the parameters are initialized, the behavior will not be identical to the vanilla model --- we will fix this below, but for now, we focus on model construction.
#
# We see that the model is a subclass of `nn.Module`, which is a subclass of Python's dataclass. The child class (written by the user) must define a `model.call(inputs)` method, that applies the function to the input, and a `model.setup()` method, that creates the modules inside this model.
#
# The module (parent) class defines two main methods: `model.apply(variables, input`, that applies the function to the input (and variables) to generate an output; and `model.init(key, input)`, that initializes the variables and returns them as a "frozen dictionary". This dictionary can contain multiple *kinds* of variables. In the example below, the only kind are parameters, which are immutable variables (that will usually get updated in an external optimization loop, as we show later). The parameters are  automatically named after the corresponding module (here, dense0, dense1, etc).  In this example, both modules are dense layers, so their parameters are a weight matrix (called 'kernel') and a bias vector.
#
# The hyper-parameters (in this case, the size of each layer) are stored as attributes of the class, and are specified when the module is constructed.

# + id="3zueDo1r0Qav"
class MLP(nn.Module):
  features: Sequence[int]
  default_attr: int = 42

  def setup(self):
    print('setup')
    self.layers = [nn.Dense(feat) for feat in self.features]

  def __call__(self, inputs):
    print('call')
    x = inputs
    for i, lyr in enumerate(self.layers):
      x = lyr(x)
      if i != len(self.layers) - 1:
        x = nn.relu(x)
    return x

# + colab={"base_uri": "https://localhost:8080/"} id="OoYDn8lX7_ZH" outputId="9aed2723-3248-417b-9a25-408ef49763d7"
key = random.PRNGKey(0)
D = 3
N = 2
x = random.normal(key, (N,D,))
layer_sizes = [3,1] # 1 hidden layer of size 3, 1 scalar output

print('calling constructor')
model = MLP(layer_sizes) # just initialize attributes of the object
print('OUTPUT')
print(model)

print('\ncalling init')
variables = model.init(key, x)  # calls setup then __call___
print('OUTPUT')
print(variables)

print('Calling apply')
y = model.apply(variables, x) # calls setup then __call___
print(y)

# + [markdown] id="5lwM1j1WksDG"
# # Compact modules
#
# To reduce the amount of boiler plate code, flax makes it possible to define a module just by writing the `call` method, avoiding the need to write a `setup` function. The corresponding layers will be created when the `init` funciton is called, so the input shape can be inferred lazily (when passed an input).

# + colab={"base_uri": "https://localhost:8080/"} id="Akq_iXXdktwb" outputId="2827db63-b8a0-4500-dd78-bd6b49f72065"
class MLP(nn.Module):
  features: Sequence[int]

  @nn.compact
  def __call__(self, inputs):
    x = inputs
    for i, feat in enumerate(self.features):
      x = nn.Dense(feat)(x)
      if i != len(self.features) - 1:
        x = nn.relu(x)
    return x

model = MLP(layer_sizes)
print(model)

params = model.init(key, x)
print(params)

y = model.apply(params, x)
print(y)

# + [markdown] id="dNiuZ54yB7Gj"
# # Explicit parameter initialization
#
# We can control the initialization of the random parameters in each submodule by specifying an init function. Below we show how to initialize our MLP to match the vanilla JAX model. We then check both methods give the same outputs.

# + id="5W_lEFsU4t04"
def make_const_init(x):
  def init_params(key, shape, dtype=jnp.float32):
    return x
  return init_params

class MLP_init(nn.Module):
  features: Sequence[int]
  params_init: Dict

  def setup(self):
    nlayers = len(self.features)
    layers = []
    for i in range(nlayers):
      W = self.params_init[f'W{i}'];
      b = self.params_init[f'b{i}'];
      weights_init = make_const_init(W)
      bias_init = make_const_init(b)
      layer = nn.Dense(self.features[i], kernel_init=weights_init, bias_init=bias_init)
      layers.append(layer)
    self.layers = layers

  def __call__(self, inputs):
    x = inputs
    for i, lyr in enumerate(self.layers):
      x = lyr(x)
      if i != len(self.layers) - 1:
        x = nn.relu(x)
    return x

# + colab={"base_uri": "https://localhost:8080/"} id="R27PhzrLY_zJ" outputId="adaba9d1-90a4-444b-95be-cef4da0768fe"
params_init = params0
model = MLP_init(layer_sizes, params_init)
print(model)

variables = model.init(key, x)
params = variables['params']
print(params)

W0 = params0['W0']
W = params['layers_0']['kernel']
assert np.allclose(W, W0)

y = model.apply(variables, x)
print(y)
assert np.allclose(y, y0)

# + [markdown] id="Rf8avaA_nGJ1"
# # Creating your own modules
#
# Now we illustrate how to create a module with its own parameters, instead of relying on composing built-in primitives. As an example, we write our own dense layer class.

# + colab={"base_uri": "https://localhost:8080/"} id="WUJ98XpSnS8F" outputId="ccdc09ef-f87f-4234-d64e-0bd583406851"
class SimpleDense(nn.Module):
  features: int # num output features for this layer
  kernel_init: Callable = nn.initializers.lecun_normal()
  bias_init: Callable = nn.initializers.zeros

  @nn.compact
  def __call__(self, inputs):
    features_in = inputs.shape[-1] # infer shape from input
    features_out = self.features
    kernel = self.param('kernel', self.kernel_init, (features_in, features_out))
    bias = self.param('bias', self.bias_init, (features_out,))
    outputs = jnp.dot(inputs, kernel) + bias
    return outputs

model = SimpleDense(features=3)
print(model)

vars = model.init(key, x)
print(vars)

y = model.apply(vars, x)
print(y)

# + [markdown] id="BJpBh933-GTW"
# # Stochastic layers
#
# Some layers may need a source of randomness. If so, we must pass them a PRNG in the `init` and `apply` functions, in addition to the PRNG used for parameter initialization. We illustrate this below using dropout. We construct two versions, one which is stochastic (for training), and one which is deterministic (for evaluation).

# + colab={"base_uri": "https://localhost:8080/"} id="tMSpLucO-Yfj" outputId="9726f1f3-ca25-4946-f9da-29692b1b034c"
class Block(nn.Module):
  features: int
  training: bool
  @nn.compact
  def __call__(self, inputs):
    x = nn.Dense(self.features)(inputs)
    x = nn.Dropout(rate=0.5)(x, deterministic=not self.training)
    return x

N = 1; D = 2;
x = random.uniform(key, (N,D))

model = Block(features=3, training=True)
key = random.PRNGKey(0)
variables = model.init({'params': key, 'dropout': key}, x)
#variables = model.init(key, x) # cannot share the rng
print('variables', variables)

# Apply stochastic model
for i in range(2):
  key, subkey = random.split(key)
  y = model.apply(variables, x, rngs={'dropout': subkey})
  print(f'train output {i}, ', y)

# Now make a deterministic version
eval_model = Block(features=3, training=False)
key = random.PRNGKey(0)
#variables = eval_model.init({'params': key, 'dropout': key}, x)
for i in range(2):
  key, subkey = random.split(key)
  y = eval_model.apply(variables, x, rngs={'dropout': subkey})
  print(f'eval output {i}, ', y)

# + [markdown] id="CHieB2aAumdg"
# # Stateful layers
#
# In addition to parameters, linen modules can contain other kinds of variables, which may be mutable as we illustrate below.
# Indeed, parameters are just a special case of variable.
# In particular, this line
# ```
# p = self.param('param_name', init_fn, shape, dtype)
# ```
# is a convenient shorthand for this:
# ```
# p = self.variable('params', 'param_name', lambda s, d: init_fn(self.make_rng('params'), s, d), shape, dtype).value
# ```
#

# + [markdown] id="EAQxx2Tu8xln"
# ## Example: counter

# + colab={"base_uri": "https://localhost:8080/"} id="BeGNa8zaut41" outputId="eb8923e0-ed62-46f9-f11b-80dec053e31a"
class Counter(nn.Module):
  @nn.compact
  def __call__(self):
    # variable(collection, name, init_fn, *init_args)
    counter1 = self.variable('counter', 'count1', lambda: jnp.zeros((), jnp.int32))
    counter2 = self.variable('counter', 'count2', lambda: jnp.zeros((), jnp.int32))
    is_initialized = self.has_variable('counter', 'count1')
    if is_initialized:
      counter1.value += 1
      counter2.value += 2
    return counter1.value, counter2.value

model = Counter()
print(model)

init_variables = model.init(key) # calls the `call` method
print('initialized variables:\n', init_variables)
counter = init_variables['counter']['count1']
print('counter 1 value', counter)

y, mutated_variables = model.apply(init_variables, mutable=['counter'])
print('mutated variables:\n', mutated_variables)
print('output:\n', y)

# + [markdown] id="1IaC2RT1v65t"
# ## Combining mutable variables and immutable parameters
#
# We can combine mutable variables with immutable parameters.
# As an example, consider a simplified version of batch normalization, which
#  computes the running mean of its inputs, and adds an optimzable offset (bias) term.
#
#

# + id="NXP19telv_Y_"
class BiasAdderWithRunningMean(nn.Module):
  decay: float = 0.99

  @nn.compact
  def __call__(self, x):
    is_initialized = self.has_variable('params', 'bias')

    # variable(collection, name, init_fn, *init_args)
    ra_mean = self.variable('batch_stats', 'mean', lambda s: jnp.zeros(s), x.shape[1:])

    dummy_mutable = self.variable('mutables', 'dummy', lambda s: 42, 0)

    # param(name, init_fn, *init_args)
    bias = self.param('bias', lambda rng, shape: jnp.ones(shape), x.shape[1:])

    if is_initialized:
      ra_mean.value = self.decay * ra_mean.value + (1.0 - self.decay) * jnp.mean(x, axis=0, keepdims=True)

    return x - ra_mean.value + bias

# + [markdown] id="x_WsMGY8xA_x"
#
# The intial variables are:
# params = (bias=1), batch_stats=(mean=0)
#
# If we pass in x=ones(N,D), the  running average becomes
# $$
# 0.99*0 + (1-0.99)*1 = 0.01
# $$
# and the output becomes
# $$
# 1 - 0.01 + 1 = 1.99
# $$
#
#
#

# + colab={"base_uri": "https://localhost:8080/"} id="dvXKCE8yxiTu" outputId="4ddb1117-32a0-481d-d8bb-876010d0e821"
key = random.PRNGKey(0)
N = 2
D = 5
x = jnp.ones((N,D))
model = BiasAdderWithRunningMean()

variables = model.init(key, x)
print('initial variables:\n', variables)
nonstats, stats = variables.pop('batch_stats')
print('nonstats', nonstats)
print('stats', stats)

# + colab={"base_uri": "https://localhost:8080/"} id="Ytr2_w9U12PT" outputId="30555a51-9b09-4ef2-8222-98c9b52e4a47"
y, mutables = model.apply(variables, x, mutable=['batch_stats'])
print('output', y)
print('mutables', mutables)

# + [markdown] id="B1g2GW3f3B-Z"
# To call the function with the updated batch stats, we have to stitch together the new mutated state with the old state, as shown below.
#
#

# + colab={"base_uri": "https://localhost:8080/"} id="cpBb21A72Bdj" outputId="d5ce7521-90c4-48a0-a180-4f02be5fe5f8"

variables = unfreeze(nonstats)
print(variables)
variables['batch_stats'] = mutables['batch_stats']
variables = freeze(variables)
print(variables)

# + [markdown] id="sa7nH74Y5-Lg"
# If we pass in x=2*ones(N,D), the running average gets updated to
# $$
# 0.99 * 0.01 + (1-0.99) * 2.0 = 0.0299
# $$
# and the output becomes
# $$
# 2- 0.0299 + 1 = 2.9701
# $$

# + colab={"base_uri": "https://localhost:8080/"} id="t2dF2si51QN5" outputId="6177ee1c-41b7-40e6-d954-d0c1c85b0180"

x = 2*jnp.ones((N,D))
y, mutables = model.apply(variables, x, mutable=['batch_stats'])
print('output', y)
print('batch_stats', mutables)

assert np.allclose(y, 2.9701)
assert np.allclose(mutables['batch_stats']['mean'], 0.0299)

# + [markdown] id="cnBmgGxOoPKU"
# # Optimization
#
# Flax has several built-in (first-order) optimizers, as we illustrate below on a random linear function. (Note that we can also fit a model defined in flax using some other kind of optimizer, such as that provided by the [optax library](https://github.com/deepmind/optax).)

# + colab={"base_uri": "https://localhost:8080/"} id="OTHgj_pMra3H" outputId="d142c2bb-c725-47e6-8cba-5b55f9f16b48"
D = 5
key = jax.random.PRNGKey(0)
params = {'w': jax.random.normal(key, (D,))}
print(params)

x = jax.random.normal(key, (D,))

def loss(params):
  w = params['w']
  return jnp.dot(x, w)

loss_grad_fn = jax.value_and_grad(loss)
v, g = loss_grad_fn(params)
print(v)
print(g)

# + id="7KdmBHa8oWFY" colab={"base_uri": "https://localhost:8080/"} outputId="3695c5f5-339c-4e27-d80e-30d100ddae66"
from flax import optim
optimizer_def = optim.Momentum(learning_rate=0.1, beta=0.9)
print(optimizer_def)

optimizer = optimizer_def.create(params)
print(optimizer)

# + colab={"base_uri": "https://localhost:8080/"} id="1JpgauX_ox_w" outputId="c649c61c-93ba-41a0-a834-2c254a53a243"
for i in range(10):
  params = optimizer.target
  loss_val, grad = loss_grad_fn(params)
  optimizer = optimizer.apply_gradient(grad)
  params = optimizer.target
  print('step {}, loss {:0.3f}, params {}'.format(i, loss_val, params))

# + [markdown] id="_ITTDWT2ECxC"
# # Worked example: MLP for MNIST
#
# We demonstrate how to fit a shallow MLP to MNIST using Flax.
# We use this function:
# https://github.com/probml/pyprobml/blob/master/scripts/fit_flax.py
# To allow us to edit this file locally (in colab), and push commits back to github, we sync this colab with github. (For details see [this colab](https://colab.research.google.com/github/probml/pyprobml/blob/master/book1/intro/colab_intro.ipynb), the cell labeled "Working with github".)
#
#

# + [markdown] id="vavamofruHS_"
# ## Import code

# + colab={"base_uri": "https://localhost:8080/"} id="L3KR7bMQWQ2A" outputId="6d6f555c-4da0-46dd-975c-ce7daeb24564"
# !ls

# + colab={"base_uri": "https://localhost:8080/"} id="3OhoJQUnuE1Q" outputId="d40b0ab5-8fb3-4286-a925-51e41be70d40"
from google.colab import drive
drive.mount('/content/drive')
# !ls /content/drive/MyDrive/ssh/

# + colab={"base_uri": "https://localhost:8080/"} id="ItaIH9dyocZA" outputId="aebdafd4-9c3c-4bb5-a9a7-f72dbe21f979"
# !rm -rf probml_tools*.*
# !wget https://raw.githubusercontent.com/probml/pyprobml/master/scripts/probml_tools.py
import probml_tools as pml

# + colab={"base_uri": "https://localhost:8080/"} id="rUYY3-N4uWYh" outputId="a9e455f0-0b49-4429-ad1d-689f5a6936cc"

# !rm -rf pyprobml
pml.git_ssh("git clone https://github.com/probml/pyprobml.git")

# + id="fjUVc9yKvFFP"
# %load_ext autoreload
# %autoreload 2

# + colab={"base_uri": "https://localhost:8080/", "height": 17} id="OQ7d5hfulh7r" outputId="744eedb7-72d1-4912-fe55-bd79c3fbd3fe"
from google.colab import files
files.view('/content/pyprobml/scripts/fit_flax.py')

# + colab={"base_uri": "https://localhost:8080/"} id="gKnFd4MFu_jE" outputId="0d715cf0-8b5c-495d-a1c4-4357d3fe1038"
import pyprobml.scripts.fit_flax as ff
ff.test()

# + [markdown] id="AMjR542vpGAQ"
# Edit the file, then commit changes.

# + colab={"base_uri": "https://localhost:8080/"} id="YJXwfqz0-_XJ" outputId="01fe932f-10f2-4486-c8bd-80ef64c06e68"
# If made any local changes to fit_flax.py, save them to github
# %cd /content/pyprobml
pml.git_ssh("git add scripts; git commit -m 'push from colab'; git push")
# %cd /content

# + [markdown] id="E_xSZi3v03pC"
# ## Data

# + colab={"base_uri": "https://localhost:8080/"} id="l4uNqjBIW0we" outputId="566ff9c6-ca6f-42a7-dbf4-abf2c78f6d53"

def process_record(batch):
  image = batch['image']
  label = batch['label']
  # flatten image to vector
  shape = image.get_shape().as_list()
  D = np.prod(shape) # no batch dimension
  image = tf.reshape(image, (D,))
  # rescale to -1..+1
  image = tf.cast(image, dtype=tf.float32)
  image = ((image / 255.) - .5) * 2.
  # convert to standard names
  return {'X': image, 'y': label}

def load_mnist(split, batch_size):
  dataset, info = tfds.load("mnist", split=split, with_info=True)
  dataset = dataset.map(process_record)
  if split=="train":
    dataset = dataset.shuffle(10*batch_size, seed=0)
  dataset = dataset.batch(batch_size)
  dataset = dataset.prefetch(tf.data.experimental.AUTOTUNE)
  dataset = dataset.cache()
  dataset = dataset.repeat()
  dataset = tfds.as_numpy(dataset) # leave TF behind
  num_examples = info.splits[split].num_examples
  return iter(dataset), num_examples

batch_size = 100
train_iter, num_train = load_mnist("train", batch_size)
test_iter, num_test = load_mnist("test", batch_size)

num_epochs = 3
num_steps = num_train // batch_size
print(f'{num_epochs} epochs with batch size {batch_size} will take {num_steps} steps')

batch = next(train_iter)
print(batch['X'].shape)
print(batch['y'].shape)

# + [markdown] id="rLiWUSjR05BQ"
# ## Model
#
#

# + id="cLwAwqd4Nzvy"
class Model(nn.Module):
  nhidden: int
  nclasses: int

  @nn.compact
  def __call__(self, x):
    if self.nhidden > 0:
      x = nn.Dense(self.nhidden)(x)
      x = nn.relu(x)
    x = nn.Dense(self.nclasses)(x) # logits
    x = nn.log_softmax(x) # log probabilities
    return x

# + [markdown] id="9JsVFGfU628j"
# ## Training loop
#

# + colab={"base_uri": "https://localhost:8080/", "height": 497} id="KDAJthPTvxI7" outputId="3f277974-c0db-4c0e-c39d-04648ae16f8f"

model = Model(nhidden = 128, nclasses=10)
rng = jax.random.PRNGKey(0)
num_steps = 200

params, history =  ff.fit_model(
    model, rng, num_steps, train_iter, test_iter, print_every=20)

display(history)

# + colab={"base_uri": "https://localhost:8080/", "height": 278} id="RWv2Sspl8EAN" outputId="88ac8b9a-9bc6-4921-8e0c-71ec0d61210d"
plt.figure()
plt.plot(history['step'], history['test_accuracy'], 'o-', label='test accuracy')
plt.xlabel('num. minibatches')
plt.legend()
plt.show()

# + id="oWe69Z51Q3Kz"

max_stars_repo_name: arpitvaghela/probml-notebooks

max_issues_repo_name: arpitvaghela/probml-notebooks

max_forks_repo_name: arpitvaghela/probml-notebooks

reason: import numpy,import jax,from jax


text:
# import necessary packages

import cv2
import sys
import numpy as np
import mediapipe as mp
import tensorflow as tf
from tensorflow.keras.models import load_model
import serial
import time
import os
import torch
sys.path.insert(0, os.path.join(os.path.dirname(__file__), '../'))
import ChessboardCalibration_master.tools.predict as predict
import object_detection

def yolo(image):
    model = torch.hub.load('ultralytics/yolov5', 'custom', path='./model/super-best.pt')
    model.conf = 0.6
    model.iou = 0.6
    img1 = image[..., ::-1]
    results = model(img1)
    results.print()
    # results.show()
    # results.save()
    results.xywh[0]
    a = results.pandas().xywh[0]
    # print(a)
    results.show()
    # print(a['confidence'][0])
    boxes = []
    return a

def get_corr(img):
    # img = cv2.imread(img_path)
    corr_data = yolo(img)
    red_corr = corr_data[corr_data['name'] == 'red']
    xcenter_pixel = red_corr['xcenter'].values[0]
    ycenter_pixel = red_corr['ycenter'].values[0]
    #67pixel = 3cm
    # time.sleep(2)
    xcenter_mm, ycenter_mm = predict.pred(xcenter_pixel, ycenter_pixel, "./robot_cam/image.jpg")
    # xcenter_mm = xcenter_pixel/62*30
    xcenter_mm = round(xcenter_mm,2)*10
    # ycenter_mm = ycenter_pixel/62*30
    ycenter_mm = round(ycenter_mm, 2)*10
    # err_x = abs(int(x)-xcenter_mm)
    # err_x = round(err_x,2)
    # err_y = abs(int(y)-ycenter_mm)
    # err_y = round(err_y,2)
    import math
    # err = math.sqrt(err_x*err_x+err_y*err_y)
    # err = round(err, 2)
    return "{} {}".format(xcenter_mm, ycenter_mm)

# ser = serial.Serial("/dev/ttyUSB0", 9600, timeout=1)

def write_data(string):
    ser.write(string.encode())
    ser.close()

# initialize mediapipe
mpHands = mp.solutions.hands
hands = mpHands.Hands(max_num_hands=1, min_detection_confidence=0.7)
mpDraw = mp.solutions.drawing_utils

# Load the gesture recognizer model
model = load_model('./gesture/mp_hand_gesture')

# Load class names
f = open('./gesture/gesture.names', 'r')
classNames = f.read().split('\n')
f.close()
print(classNames)

# Initialize the webcam
cap1 = cv2.VideoCapture(1)
cap2 = cv2.VideoCapture(0)

while True:
    # Read each frame from the webcam
    ret1, frame = cap1.read()
    ret2, img = cap2.read()
    if ret1:
        x, y, c = frame.shape

        # Flip the frame vertically
        frame = cv2.flip(frame, 1)
        framergb = cv2.cvtColor(frame, cv2.COLOR_BGR2RGB)

        # Get hand landmark prediction
        result = hands.process(framergb)

        # print(result)

        className = ''

        # post process the result
        if result.multi_hand_landmarks:
            landmarks = []
            for handslms in result.multi_hand_landmarks:
                for lm in handslms.landmark:
                    # print(id, lm)
                    lmx = int(lm.x * x)
                    lmy = int(lm.y * y)

                    landmarks.append([lmx, lmy])

                # Drawing landmarks on frames
                mpDraw.draw_landmarks(frame, handslms, mpHands.HAND_CONNECTIONS)

                # Predict gesture
                prediction = model.predict([landmarks])
                # print(prediction)
                classID = np.argmax(prediction)
                className = classNames[classID]

        # show the prediction on the frame
        cv2.putText(frame, className, (10, 50), cv2.FONT_HERSHEY_SIMPLEX,
                            1, (0,0,255), 2, cv2.LINE_AA)
        if cv2.waitKey(33) == ord("a"):
            if (className == "004"):
                print("hi")
                # cv2.imwrite("./robot_cam/image.jpg", cv2.imread('./a.jpg'))
                cv2.imwrite("./robot_cam/image.jpg", img)

                str = get_corr(img)
                print(str)
            # arduino.write_data(className)
            else:
                print("gui di: ", className)
                write_data(className)

        # Show the final output
        cv2.imshow("human_cam", frame)
    if ret2:
        cv2.imshow("object_cam", img)
    if cv2.waitKey(1) == ord('q'):
        break

# release the webcam and destroy all active windows
cap1.release()
cap2.release()

cv2.destroyAllWindows()


text:
import os
os.add_dll_directory("C:\\Program Files\\Thorlabs\\Kinesis")

from ctypes import c_int, c_double, byref, pointer
import numpy as np
import time
import keyboard

from pyrolab.drivers.motion.z825b import Z825B
from pyrolab.drivers.motion.prm1z8 import PRM1Z8

help_txt = """\nControls\n--------\n
    \tleft arrow - jog move left
    \tright arrow - jog move right
    \tdown arrow - jog move down
    \tup arrow - jog move up
    \tc - jog clockwise rotation
    \tx - jog counterclockwise rotation
    \ta - continuous move left
    \td - continuous move right
    \ts - continuous move down
    \tw - continuous move up
    \tshift + j - jog step size (linear motors)
    \tshift + g - jog step size (rotational motor)
    \tshift + k - set velocity (mm/s)
    \tspacebar - emergency stop all
    \to - home linear motors
    \th - help
    \tq - quit
    """

x_mot = Z825B("27504851", home=False)
y_mot = Z825B("27003497", home=False)
r_mot = PRM1Z8("27003366", home=False)
motors = [x_mot, y_mot, r_mot]

conversion_matrix = np.array([[1]])

point1 = [0,0]
point2 = [0,0]
point3 = [0,0]

x_mot_moving = False
y_mot_moving = False
r_mot_moving = False

def stop_cont_jog(motor, move_type):
    if move_type == "continuous":
        motor.stop()
    if motor == x_mot:
        global x_mot_moving
        x_mot_moving = False
    if motor == y_mot:
        global y_mot_moving
        y_mot_moving = False
    if motor == r_mot:
        global r_mot_moving
        r_mot_moving = False

def set_velocity():
    velocity = float(input('New velocity (device units):'))
    for m in motors:
        m.velocity = velocity

def set_jog_step_linear_input():
    step = float(input('New Jog Step (mm):'))
    x_mot.jog_step_size = step
    y_mot.jog_step_size = step

def set_jog_step_linear(step_size):
    x_mot.jog_step_size = step_size
    y_mot.jog_step_size = step_size

def set_jog_step_rotational():
    step = float(input('New Jog Step (degrees):'))
    r_mot.jog_step_size = step

def stop_all():
    for m in motors:
        m.stop()

def move_step(motor, direction):
    if motor == x_mot:
        global x_mot_moving
        if x_mot_moving == False:
            motor.jog(direction)
            #x_mot_moving = True
    if motor == y_mot:
        global y_mot_moving
        if y_mot_moving == False:
            motor.jog(direction)
            #y_mot_moving = True
    if motor == r_mot:
        global r_mot_moving
        if r_mot_moving == False:
            motor.jog(direction)
            #r_mot_moving = True

def move_cont(motor, direction):
    if motor == x_mot:
        global x_mot_moving
        if x_mot_moving == False:
            motor.move_continuous(direction)
            x_mot_moving = True
    if motor == y_mot:
        global y_mot_moving
        if y_mot_moving == False:
            motor.move_continuous(direction)
            y_mot_moving = True
    if motor == r_mot:
        global r_mot_moving
        if r_mot_moving == False:
            motor.move_continuous(direction)
            r_mot_moving = True

def help_me():
    global help_txt
    print(help_txt)

def set_point_1():
    global point1, dataScanner
    point1 = [x_mot.get_position(), y_mot.get_position()]
    print("Point 1 set to (" + str(point1[0]) + "," + str(point1[1]) + ")")

def set_point_2():
    global point2
    point2 = [x_mot.get_position(), y_mot.get_position()]
    print("Point 2 set to (" + str(point2[0]) + "," + str(point2[1]) + ")")

def set_point_3():
    global point3
    point3 = [x_mot.get_position(), y_mot.get_position()]
    print("Point 3 set to (" + str(point3[0]) + "," + str(point3[1]) + ")")

def keyloop_for_calibration():
    keyboard.on_press_key('left arrow', lambda _:move_step(x_mot, "backward"))
    keyboard.on_release_key('left arrow', lambda _:stop_cont_jog(x_mot, "step"))
    keyboard.on_press_key('right arrow', lambda _:move_step(x_mot, "forward"))
    keyboard.on_release_key('right arrow', lambda _:stop_cont_jog(x_mot, "step"))
    keyboard.on_press_key('up arrow', lambda _:move_step(y_mot, "forward"))
    keyboard.on_release_key('up arrow', lambda _:stop_cont_jog(y_mot, "step"))
    keyboard.on_press_key('down arrow', lambda _:move_step(y_mot, "backward"))
    keyboard.on_release_key('down arrow', lambda _:stop_cont_jog(y_mot, "step"))
    keyboard.on_press_key('a', lambda _:move_cont(x_mot, "backward"))
    keyboard.on_release_key('a', lambda _:stop_cont_jog(x_mot, "continuous"))
    keyboard.on_press_key('d', lambda _:move_cont(x_mot, "forward"))
    keyboard.on_release_key('d', lambda _:stop_cont_jog(x_mot, "continuous"))
    keyboard.on_press_key('w', lambda _:move_cont(y_mot, "forward"))
    keyboard.on_release_key('w', lambda _:stop_cont_jog(y_mot, "continuous"))
    keyboard.on_press_key('s', lambda _:move_cont(y_mot, "backward"))
    keyboard.on_release_key('s', lambda _:stop_cont_jog(y_mot, "continuous"))
    keyboard.add_hotkey('shift + j', set_jog_step_linear_input)
    keyboard.add_hotkey('shift + g', set_jog_step_rotational)
    keyboard.add_hotkey('shift + k', set_velocity)
    keyboard.on_press_key('space', lambda _:stop_all())
    keyboard.on_press_key('h', lambda _:help_me())
    keyboard.on_press_key('o', lambda _:home_motors())
    keyboard.add_hotkey('shift + 1', set_point_1)
    keyboard.add_hotkey('shift + 2', set_point_2)
    keyboard.add_hotkey('shift + 3', set_point_3)

    print(help_txt)

    while True:
        if keyboard.is_pressed('q'):
            print("QUIT")
            break

def keyloop(quitKey = 'q'):

    keyboard.on_press_key('left arrow', lambda _:move_step(x_mot, "backward"))
    keyboard.on_release_key('left arrow', lambda _:stop_cont_jog(x_mot, "step"))
    keyboard.on_press_key('right arrow', lambda _:move_step(x_mot, "forward"))
    keyboard.on_release_key('right arrow', lambda _:stop_cont_jog(x_mot, "step"))
    keyboard.on_press_key('up arrow', lambda _:move_step(y_mot, "forward"))
    keyboard.on_release_key('up arrow', lambda _:stop_cont_jog(y_mot, "step"))
    keyboard.on_press_key('down arrow', lambda _:move_step(y_mot, "backward"))
    keyboard.on_release_key('down arrow', lambda _:stop_cont_jog(y_mot, "step"))
    keyboard.on_press_key('c', lambda _:move_step(r_mot, "forward"))
    keyboard.on_release_key('c', lambda _:stop_cont_jog(r_mot, "step"))
    keyboard.on_press_key('x', lambda _:move_step(r_mot, "backward"))
    keyboard.on_release_key('x', lambda _:stop_cont_jog(r_mot, "step"))
    keyboard.on_press_key('a', lambda _:move_cont(x_mot, "backward"))
    keyboard.on_release_key('a', lambda _:stop_cont_jog(x_mot, "continuous"))
    keyboard.on_press_key('d', lambda _:move_cont(x_mot, "forward"))
    keyboard.on_release_key('d', lambda _:stop_cont_jog(x_mot, "continuous"))
    keyboard.on_press_key('w', lambda _:move_cont(y_mot, "forward"))
    keyboard.on_release_key('w', lambda _:stop_cont_jog(y_mot, "continuous"))
    keyboard.on_press_key('s', lambda _:move_cont(y_mot, "backward"))
    keyboard.on_release_key('s', lambda _:stop_cont_jog(y_mot, "continuous"))
    keyboard.add_hotkey('shift + j', set_jog_step_linear_input)
    keyboard.add_hotkey('shift + g', set_jog_step_rotational)
    keyboard.add_hotkey('shift + k', set_velocity)
    keyboard.on_press_key('space', lambda _:stop_all())
    keyboard.on_press_key('h', lambda _:help_me())
    keyboard.on_press_key('o', lambda _:home_motors())
    #keyboard.on_press_key('d', get_data)

    print(help_txt)

    while True:
        if keyboard.is_pressed(quitKey):
            print("QUIT")
            break

def home_motors():
    print("Homing motors...")
    x_mot.go_home()
    y_mot.go_home()
    print("Done homing")

def go_to_GDS_Coordinates(x, y, conversion_matrix):
    oldPoint = np.array([[x], [y], [1]])
    newPoint = conversion_matrix @ oldPoint
    x_mot.move_to(newPoint[0][0])
    y_mot.move_to(newPoint[1][0])

def go_to_chip_coordinates(x, y):
    x_mot.move_to(x)
    y_mot.move_to(y)

def go_to_chip_coordinates_y(y):
    y_mot.move_to(y)

def go_to_chip_coordinates_x(x):
    x_mot.move_to(x)

def go_to_circuit(circuit):
    pos = circuit.location
    go_to_GDS_Coordinates(float(pos[0]), float(pos[1]))

def calibrate(circuits):
    global point1, point2, point3, conversion_matrix
    home_motors()
    print("Starting calibration process...")
    print("Move " + circuits[0].name + " into Crosshairs, then press shift + 1")
    print("Then Move " + circuits[1].name + " into Crosshairs, then press shift + 2")
    print("Then Move " + circuits[2].name + " into Crosshairs, then press shift + 3")
    keyloop_for_calibration()
    location1 = circuits[0].location
    location2 = circuits[1].location
    location3 = circuits[2].location
    print("GDS locations:")
    print(location1)
    print(location2)
    print(location3)
    print("Chip locations:")
    print(point1)
    print(point2)
    print(point3)

    GDS = np.array([[float(location1[0]), float(location1[1]), 1, 0, 0, 0], [0, 0, 0, float(location1[0]), float(location1[1]), 1], [float(location2[0]), float(location2[1]), 1, 0, 0, 0], [0, 0, 0, float(location2[0]), float(location2[1]), 1], [float(location3[0]), float(location3[1]), 1, 0, 0, 0], [0, 0, 0, float(location3[0]), float(location3[1]), 1]])

    CHIP = np.array([[point1[0]],[point1[1]],[point2[0]],[point2[1]],[point3[0]],[point3[1]]])

    a = np.linalg.inv(GDS) @ CHIP

    conversion_matrix = np.array([[a[0][0], a[1][0], a[2][0]], [a[3][0], a[4][0], a[5][0]], [0, 0, 1]])

def move_right():
    move_step(x_mot, 'forward')

def move_left():
    move_step(x_mot, 'backward')

def move_up():
    move_step(y_mot, 'forward')

def move_down():
    move_step(y_mot, 'backward')

def get_y_position():
    return y_mot.get_position()

def get_x_position():
    return x_mot.get_position()

max_stars_repo_name: BYUCamachoLab/autogator

max_issues_repo_name: BYUCamachoLab/autogator

max_forks_repo_name: BYUCamachoLab/autogator


text:
from matplotlib import pyplot as plt
from os import path
import numpy as np
import time

import sys
sys.path.append(path.split(path.dirname(__file__))[0])

from plato.backend.units import Units

from plato.backend.datasources.branch_training_trace.datasource import BranchTrainingDatasource
from plato.backend.processors.branch_training_heatmap.adapter import BranchTrainingHeatMapAdapter
from plato.backend.processors.branch_training_heatmap.generator import BranchTrainingHeatMapGenerator

if len(sys.argv) == 1:
    filename = path.join(path.dirname(__file__),'test-branch-training-trace.hdf5')
else:
    filename = sys.argv[1]

def plot_heatmap(hm):
    plt.figure()
    plt.imshow(hm, cmap='hot', interpolation='nearest', aspect='auto')

# Load source data
branch_hm_data = BranchTrainingDatasource(filename)
print('\nstats', branch_hm_data.stats)
bphma = BranchTrainingHeatMapAdapter(branch_hm_data)
print('\nnum events', bphma.num_events)

bin_size = 200000

bphmg = BranchTrainingHeatMapGenerator(bphma, ['thrash_1', 'd_weight'], bin_size)

hm, table_means, row_means = bphmg.generate_2d_heatmap_with_profiles(0, len(branch_hm_data.ddf_branch_training_events)-1, Units.BRANCH_TRAINING_INDEX, 'd_weight')

print('max {} at {}'.format(hm.max(), hm.argmax()))
print('min {} at {}'.format(hm.min(), hm.argmin()))

hm, table_means, row_means = bphmg.generate_2d_heatmap_with_profiles(0, len(branch_hm_data.ddf_branch_training_events)-1, Units.BRANCH_TRAINING_INDEX, 'd_weight', allow_bins=False)

print('max {} at {}'.format(hm.max(), hm.argmax()))
print('min {} at {}'.format(hm.min(), hm.argmin()))

plot_heatmap(hm)
plt.title('d_weight')

# Show all the images (blocking)
sys.stdout.flush()

import os
plt.show(block=bool(int(os.environ.get('BLOCKING_SHOW', 1))))


text:
"""M3C 2017 Homework 1 - Lawrence Stewart CID 00948972 -email ls3914@ic.ac.uk
"""
import numpy as np
import matplotlib.pyplot as plt
import networkx as nx

def rw2d(Nt,M,a=0,b=0):
    """Input variables
    Nt: number of time steps
    M: Number of simulations
    a,b: bias parameters
    """

    #Ensure the constraints on a,b given in the question
    assert a>=-1, "a must be greater than or equal to -1"
    assert 1>=b, "b must be less than or equal to 1"

    # Generate Nt randomly selected values for both the X and Y co-ordinate with x in {-1,1+a} and y in {-1,1-b}
    Xmat=np.random.choice(np.array([-1,1+a]),[Nt+1,M])
    Ymat=np.random.choice(np.array([-1,1-b]),[Nt+1,M])

    #Set initial conditions for Nt=0:
    Xmat[0,:]=0.
    Ymat[0,:]=0.

    # If errors occur change back to this:
    # Xmat2=np.cumsum(Xmat,axis=0)
    # Ymat2=np.cumsum(Ymat,axis=0)

    Xmat=np.cumsum(Xmat,axis=0)
    Ymat=np.cumsum(Ymat,axis=0)

    #Generate <X>,<Y>:
    X_=np.mean(Xmat,axis=1)
    Y_=np.mean(Ymat,axis=1)

    #Generate <X^2>, <Y^2>
    X2_=np.mean(np.multiply(Xmat,Xmat),axis=1)
    Y2_=np.mean(np.multiply(Ymat,Ymat),axis=1)

    #Generate <XY>
    XY_=np.mean(np.multiply(Xmat,Ymat),axis=1)

    return X_,Y_, X2_, Y2_, XY_

def rwnet1(H,Hf,a=0,display=False):
    """Input variables
    H: Height at which new nodes are initially introduced
    Hf: Final network height
    a: horizontal bias parameter, should be -1, 0, or 1
    display: figure displaying the network is created when true
    Output variables
    X,Y: Final node coordinates
    output: a tuple containing any other information you would
    like the function to return, may be left empty
    """

    #Principlal idea behind method: At each random walk, calculate the number of steps in any direction that can be done at once
    # without violating any of the conditions. Vectorise the process for all of these steps, then re-evaluate.

    #Let output be the y co-ordinate if the highest node in the network
    output=np.array([0])

    #Ensure the constraints on a, given in the question
    assert a in [-1,0,1], "Input error - a must be either 1, 0, or -1"

    #Ensure constraints placed on H, given in the question
    assert H>=0, "Input error - H must be positive"

    #Ensure the constraints placed on Hf
    assert Hf>1, "Hf must be an integer greater than 1"
    assert H+3>=Hf, "Hf must be less than or equal to H-3"

    #Variable to check that no value has been placed above Hf
    place_check=False

    #Variable to check that node has not violated either of the conditions for continued random walk
    walk_check=False

    #List of the node co-ordinates (x,y), where we start with the origin node (0,0)
    X,Y=np.array([0]),np.array([0])

    #Calculate the sqrt(1^2+(1+alpha)^2) which will the upper bound for a single movement step in optimal direction during a walk
    maxstep=np.sqrt(1+(1+a)**2)

    while place_check==False:

        #Initiate new node placement, x y are the current co-ordinates of the random walk:
        x_cur,y_cur=0,H

        #Reset Walk_check to False, for we are starting a new walk:
        walk_check=False

        #If walk_check remains false begin a new random walk:
        while walk_check==False:

            #Calculate the distance to the closest neighbour and to the y axis
            xdist,ydist=X-x_cur,Y-y_cur
            node_dist=np.sqrt(np.multiply(xdist,xdist)+np.multiply(ydist,ydist)) #This generates the distance of current point from all nodes
            min_ind=np.argmin(np.array([min(node_dist),y_cur]))

            #If min_ind=0 then the closest thing is a node and if min_ind=1 then the closest thing is the y axis

            #set the distance (if a node is the closest)
            if min_ind==0:
                d=min(node_dist)

            #set distance (if the y=0 axis is the closest)
            else:
                d=y_cur

            #Calculate the maximum steps we can do as one process safely i.e floor( d/ maxstep) [note int() rounds down]
            nskip=int(d/maxstep)

            #Check if nskip=0 (in which case we are have finished the walk)
            if nskip==0:

                #Add node co-ordinates
                X=np.append(X,x_cur)
                Y=np.append(Y,y_cur)

                #complete the walk:
                walk_check=True

                #Update the height of the network:
                output=np.append(output,max(Y))

                #If termination condition has been reached - end simulation
                if y_cur>=Hf:
                    place_check=True

            #Perform random walk for nskip iterations:

            #1) As we are doing nskip steps in one we can sample from binomial distribution for both X and Y and
            #sum the total successes, where p success =1/2  (a success will be the # of -1's while a fail will be the number of 1+a or 1-b's)

            x_rand=sum(np.random.binomial(1,0.5,nskip))
            y_rand=sum(np.random.binomial(1,0.5,nskip))

            #we have x_rand numbers of -1 steps and (nskip-xrand) number of (1+a) steps
            x_cur=x_cur+((1+a)*(nskip-x_rand)-x_rand)

            #we have y_rand number of -1 steps and (n_skip-y_rand) number of (1-b=0) steps
            y_cur=y_cur-y_rand

    #Create Plots if User Desires
    if display==True:
        plt.figure()
        plt.suptitle('Lawrence Stewart - Created Using rwnet1().')
        plt.scatter(X,Y,c=X,cmap='viridis_r')
        plt.xlabel("x co-ordinates")
        plt.ylabel("y co-ordinates")
        plt.title("Plot of Generated Network with a=%i"%a )
        plt.grid()
        plt.show()

    return X,Y,output

def rwnet2(L,H,Hf,a=0,display=False):
    """Input variables
    L: Walls are placed at X = +/- L
    H: Height at which new nodes are initially introduced
    Hf: Final network height
    a: horizontal bias parameter, should be -1, 0, or 1
    display: figure displaying the network is created when true
    Output variables
    X,Y: Final node coordinates
    output: a tuple containing any other information you would
    like the function to return, may be left empty
    """
    #Principlal idea behind method: At each random walk, calculate the number of steps in any direction that can be done at once
    # without violating any of the conditions. Vectorise the process for all of these steps, then re-evaluate.

    #Ensure the constraints on a, given in the question
    assert a in [-1,0,1], "Input error - a must be either 1, 0, or -1"

    #Ensure constraints placed on H, given in the question
    assert H>=0, "Input error - H must be positive"

    #Ensure the constraints placed on Hf
    assert Hf>1, "Hf must be an integer greater than 1"
    assert H+3>=Hf, "Hf must be less than or equal to H-3"

    #Let output be the y co-ordinate if the highest node in the network
    output=np.array([0])

    #Variable to check that no value has been placed above Hf
    place_check=False

    #Variable to check that node has not violated either of the conditions for continued random walk
    walk_check=False

    #List of the node co-ordinates (x,y), where we start with the origin node (0,0)
    X,Y=np.array([0]),np.array([0])

    #Calculate the sqrt(1^2+(1+alpha)^2) which will the upper bound for a single movement step in optimal direction during a walk
    maxstep=np.sqrt(1+(1+a)**2)

    while place_check==False:

        #Initiate new node placement, x y are the current co-ordinates of the random walk:
        x_cur,y_cur=0,H

        #Reset Walk_check to False, for we are starting a new walk:
        walk_check=False

        #If walk_check remains false begin a new random walk:
        while walk_check==False:

            #Calculate the distance to the closest neighbour and to the y axis
            xdist,ydist=X-x_cur,Y-y_cur
            node_dist=np.sqrt(np.multiply(xdist,xdist)+np.multiply(ydist,ydist)) #This generates the distance of current point from all nodes
            min_ind=np.argmin(np.array([min(node_dist),y_cur]))

            #If min_ind=0 then the closest thing is a node and if min_ind=1 then the closest thing is the y axis

            #set the distance (if a node is the closest)
            if min_ind==0:
                d=min(node_dist)

            #set distance (if the y=0 axis is the closest)
            else:
                d=y_cur

            #calculate the distance to the walls
            walld=min(abs(L-x_cur),abs(-L-x_cur))

            #Calculate the minimum number of steps that we have until hitting a wall
            wallskip=int(walld/max( 1,abs(1+a) ) )

            #take into account walls:

            #Calculate the maximum steps we can do without hitting Y=0 or another node
            nskip=int(d/maxstep)

            #Check if nskip=0 (in which case we are have finished the walk)
            if nskip==0:

                #Add node co-ordinates
                X=np.append(X,x_cur)
                Y=np.append(Y,y_cur)

                #complete the walk:
                walk_check=True

                #Update the height of the network:
                output=np.append(output,max(Y))

                #If full termination condition has been reached - end simulation
                if y_cur>=Hf:
                    place_check=True
                break

            #skip is the number of steps we can do at once
            skip=max(1,min(nskip,wallskip))

            #Perform random walk for skip iterations:

            #1) As we are doing skip steps in one we can sample from binomial distribution for both X and Y and
            #sum the total successes, where p success =1/2  (a success will be the # of -1's while a fail will be the number of 1+a or 1-b's)

            x_rand=sum(np.random.binomial(1,0.5,skip))
            y_rand=sum(np.random.binomial(1,0.5,skip))

            #we have x_rand numbers of -1 steps and (skip-xrand) number of (1+a) steps
            x_cur=x_cur+((1+a)*(skip-x_rand)-x_rand)

            #we have y_rand number of -1 steps and (skip-y_rand) number of (1-b=0) steps
            y_cur=y_cur-y_rand

            #check to see if node has gone past the walls:
            if x_cur<-L:
                x_cur=-L

            if x_cur>L:
                x_cur=L

    #Create Plots if User Desires
    if display==True:
        plt.figure()
        plt.suptitle('Lawrence Stewart - Created Using rwnet2().')
        plt.scatter(X,Y,c=X,cmap='viridis_r')
        plt.xlabel("x co-ordinates")
        plt.ylabel("y co-ordinates")
        plt.title("Plot of Generated Network with a=%i with walls at L=%i"%(a,L) ) #This line could need work
        plt.grid()
        plt.show()

    return X,Y,output

def analyze(H,Hf):
    """ Analysis:

    hw11.png, (which should contain a subplot of two graphs) shows the rate at which the Network Height Increases, with one plot for
    the unbiased scenario a=0, and the other plot for the biased scenario a=1. The independant variable is the current height of the network/number of
    nodes placed. Both plots show this for walls setup at L=30, 150 and infinity. Examining the first graph (a=0 unbiased), it can be seen that location of the walls
    has little influence on how the rate of increase of height of the network. It is also worth noting that for each of the 3 values of L, we have that the network
    stops being built after approximately 2500 iterations (where an iteration is the placement of a node.)

    However, when we examine the biased a=1 example, we see very different behaviour for different values of L. For L=30 we see that the network reaches height Hf
    within a very fast time -approximately 150 iterations. For L=150 and L=inf, we similary see the network is reaches termination height much faster, after
    approximately 900 and 1500 iterations, respectively (compare these three values to the 2500 iterations required to reach Hf in the unbiased case).

    Apart from an initial spike (in which the first few nodes are placed), the rate of growth for the unbiased test (a=0) is approx 0.07 for each value of L.
    For the biased case (a=1) we see that for L=30, the rate is approximately 1.3 on average,  implying that the network height growing very fast. For L=150, the
    rate is slightly lower at roughly 1.1, but then steadily decays with iterations. For L=infinity, we do not see the same high initial rate as with L=30,150,
    and the rate remains fairly constant (however it still takes far less iterations to reach Hf).

    This should be expected, for the expected position for a unbiased 1D random walk is the starting position. In our system, the starting point is (0,H).
    Considering only the movement in the X direction, with a=0 we should expect the majority of nodes to be near X=0. This means it is far more unlikely that nodes
    will reach X=L. Hence, the different distances for wall placements do not have much effect on the rate of increase of height of the network. This leaves the
    simulations with L=30,150,inf to have similar results for the unbiased walk, as they are undergoing a similar stochiastic process.

    With the bias added (a=1), it is expected for the majority of nodes to be shifted to the right. With L=30, it is very likely that nodes will reach the wall.
    Similary (but less likely), for L=150. The clumping of nodes close to the wall results in a higher chance of a node finishing the random walk earlier, which
    results in a higher chance of the node increasing the height of the network (as nodes can either move down or stay at the same height). With L=inf we do not
    see the same explosive growth of the network as no walls are present. However, the spread of nodes will be biased to one direction, hence the networks height will
    still grow faster than in the unbiased case. Furthermore, the distance that is required for a node to link and finish the random walk is larger with a=1 dstar=sqrt(5)
    compared to a=0 (dstar= sqrt(2)).

    hw12.png shows the height of the network vs the number of nodes placed for the 6 different combinations of a and L. As expected, for the reasons discussed above,
    we see that L=30 a=1 has the largest rate of increase of height, followed by L=150 a=1. The figure helps highlight the difference between the L=inf a=1 case vs the
    unbiased simulations. We see that despite the lack of walls, the biased network still grows at a rate substantially larger than the unbiased cases (for reasons
    discussed above.)

    """

    #Create plots for heights for all combinations for h

    print("------Generating Graphs -----",end='\r')
    _,_,Ha0Linf=rwnet1(H,Hf,0)
    print("-----Built Graph 1/6-----",end='\r')

    _,_,Ha1Linf=rwnet1(H,Hf,1)
    print("-----Built Graph 2/6-----",end='\r')

    _,_,Ha0L30=rwnet2(30,H,Hf,0)
    print("-----Built Graph 3/6-----",end='\r')

    _,_,Ha1L30=rwnet2(30,H,Hf,1)
    print("-----Built Graph 4/6-----",end='\r')

    _,_,Ha0L150=rwnet2(150,H,Hf,0)
    print("-----Built Graph 5/6-----",end='\r')

    _,_,Ha1L150=rwnet2(150,H,Hf,1)
    print("-----Built Graph 6/6-----",end='\r')

    print("-----Process Complete-----",end='\r')

    #Generate the number of nodes for each of the simulations
    nodes_Ha0Linf=np.linspace(1,np.size(Ha0Linf),np.size(Ha0Linf))
    nodes_Ha1Linf=np.linspace(1,np.size(Ha1Linf),np.size(Ha1Linf))
    nodes_Ha0L30=np.linspace(1,np.size(Ha0L30),np.size(Ha0L30))
    nodes_Ha1L30=np.linspace(1,np.size(Ha1L30),np.size(Ha1L30))
    nodes_Ha0L150=np.linspace(1,np.size(Ha0L150),np.size(Ha0L150))
    nodes_Ha1L150=np.linspace(1,np.size(Ha1L150),np.size(Ha1L150))

    #Plotting
    plt.figure(figsize=(14, 6))
    plt.suptitle('Lawrence Stewart - Created Using analyze().')

    #Plot a=0 graphs for H/#nodes i.e the rate with #of nodes
    plt.subplot(121)
    plt.plot(nodes_Ha0L30,Ha0L30/nodes_Ha0L30, 'g-', label="L = 30",alpha=0.7)
    plt.plot(nodes_Ha0L150,Ha0L150/nodes_Ha0L150, 'r-', label="L = 150",alpha=0.7)
    plt.plot(nodes_Ha0Linf,Ha0Linf/nodes_Ha0Linf, 'c-', label="L = Inf",alpha=0.7)
    plt.xlabel("Number of Nodes Placed")
    plt.ylabel("Height of Network /Number of Nodes")
    plt.title("Rate at which Network Height Increases for varied L, a=0")
    ax = plt.gca()
    ax.set_facecolor('#D9E6E8')
    plt.grid('on')
    plt.legend()

    # #plot a=1 graphs for H/#nodes i.e the rate with #of nodes
    plt.subplot(122)
    plt.plot(nodes_Ha1Linf,Ha1Linf/nodes_Ha1Linf, 'c-', label="L = Inf")
    plt.plot(nodes_Ha1L150,Ha1L150/nodes_Ha1L150, 'r-', label="L = 150")
    plt.plot(nodes_Ha1L30,Ha1L30/nodes_Ha1L30, 'g-', label="L = 30")
    plt.xlabel("Number of Nodes Placed")
    plt.ylabel("Height of Network /Number of Nodes")
    plt.title("Rate at which Network Height Increases for varied L, a=1")
    plt.grid('on')
    ax = plt.gca()
    ax.set_facecolor('#D9E6E8')
    plt.legend()
    plt.show()

    #Plot the heights vs the number of nodes for all combinations:
    plt.figure(figsize=(9, 7))
    plt.suptitle('Lawrence Stewart - Created Using analyze().')
    plt.plot(nodes_Ha0L30,Ha0L30, label="L = 30, a=0",alpha=0.7)
    plt.plot(nodes_Ha0L150,Ha0L150, label="L = 150, a=0",alpha=0.7)
    plt.plot(nodes_Ha0Linf,Ha0Linf, label="L = Inf, a=0",alpha=0.7) #need to edit the colors
    plt.plot(nodes_Ha1Linf,Ha1Linf, label="L = Inf, a=1",alpha=0.7)
    plt.plot(nodes_Ha1L150,Ha1L150, label="L = 150, a=1",alpha=0.7)
    plt.plot(nodes_Ha1L30,Ha1L30, label="L = 30, a=1",alpha=0.7)
    plt.xlabel("Number of Nodes Placed")
    plt.ylabel("Height of Network")
    plt.title("Development of Network Height with Node Placement")
    ax = plt.gca()
    ax.set_facecolor('#D9E6E8')
    plt.grid('on')
    plt.legend()
    plt.plot()
    plt.show()
    return None

def network(X,Y,dstar,display=False,degree=False):
    """ Input variables
    X,Y: Numpy arrays containing network node coordinates
    dstar2: Links are placed between nodes within a distance, d<=dstar of each other
        and dstar2 = dstar*dstar
    display: Draw graph when true
    degree: Compute, display and return degree distribution for graph when true
    Output variables:
    G: NetworkX graph corresponding to X,Y,dstar
    D: degree distribution, only returned when degree is true
    """

    """
    The function network() takes in co-ordinates (x,y) in the form of two numpy arrays X,Y. The function then loops through the arrays and compares X[i] to
    all X[j] where j>i (and does the same for Y) checking for the condition that d((x_i,y_i),(x_j,y_j))=<dstar. If this condition is met, the nodes will be linked and
     a tuple denoting the link will be added to a list that contains all the edges to be formed. From here the function uses NetworkX's inbuilt function add_edges_from()
     to build the graph. Finally plots are returned for the graph and degree distribution if the user requires.

    Analysis: It is evident that the network's maximum degree is constrained to a small value. This is because as soon as a node becomes within a proximity
    of dstar of another node, the random walk is terminated and the a connection is formed between the two nodes. Suppose for example some nodes have formed a
    "cluster" of some form. The probability that a new node penetrates into this cluster is far smaller than the probability that the node will link to
    the outside of the cluster. Inductively we can see that it is very unlikely that we will achieve high degrees in the network.

    A simple solution would be to slightly modify the function network(). Once a node is placed and a link is formed with another node
    we could then form links with all nodes in a ball of radius R (where R would be an additional input of network()). A large R value would vastly increase the maximum
    degree reached.

    A more technical change would be to alter the functionality of rwnet1 and rwnet2, so that both functions take in an additional parameter Nlink, where Nlink
    is the number of nodes that a node taking a random walk must come with a distance of dstar of to terminate the walk, i.e placing the node in the network and
    forming links. This would increase the probability of nodes penetrating deeper into clusters, and reduce the probability of the nodes attatching to the outside
    of these clusters. A problem with this method would be the placements of the first nodes, which would all end up on the Y=0 axis. A possible solution to this would
    be to introduce a dynamic value of Nlink. This could be either an array or possibly a function that maps the current iteration number to a value for Nlink. This would
    allow the simulation to initially run exactly the same is in our current code, i.e. Nlink has an initial value of 1 for the first K iterations however, after more nodes
    have been placed Nlink increases to allow nodes to form connections of higher degrees.

    """
    #Create the graph
    G=nx.Graph()

    #preallocate the edges list (and remove possible duplicates in input)
    edges=[]
    X,Y=(zip(*list(set(list(zip(X,Y) ) ))))
    X,Y=np.array(X),np.array(Y)

    for i in range(len(X)-1):

        #Calculate the distance

        tempX=X[i+1:]-X[i]
        tempY=Y[i+1:]-Y[i]
        tempd=np.sqrt(np.multiply(tempX,tempX)+np.multiply(tempY,tempY))

        #find all nodes within d*
        ind=np.where(tempd<=dstar)[0]+i+1

        #add the edges to the list
        edges_to_add=[(i,j) for j in ind] # old one list(zip(i*np.ones(np.size(ind)), ind))
        edges=edges+edges_to_add

    #Flatten list
    checker=list(sum(edges, ()))

    #Take all the nodes that have edges
    nodes_to_add=np.unique(checker)

    #add them with positions
    for i in range(np.size(X)):
        if i in checker:
            G.add_node(i, pos = (X[i], Y[i]))

    G.add_edges_from(edges)

    #Plot graph if required
    if display==True:
        plt.figure(figsize=(9, 7))
        plt.suptitle('Lawrence Stewart - Created using network().')
        # pos=nx.get_node_attributes(G,'pos')
        nx.draw_networkx(G, pos = nx.get_node_attributes(G, 'pos'), with_labels = False, node_size =20, width = 2., font_size = 6,node_color='b',prog='dot')
        # nx.draw(G,pos,node_shape='.',node_size=100,
        # prog='dot') #add this in when done  cmap=plt.cm.Blues,node_color=range(len(G)),prog='dot'
        node_color=range(len(G)),
        ax = plt.gca()
        ax.set_facecolor('#D9E6E8')
        plt.title("Graph of Network")
        plt.xlabel("x")
        plt.ylabel("y")
        plt.grid('on')
        plt.show()

    #Plot degree distribution if required
    if degree==True:
        plt.figure(figsize=(9, 7))
        plt.suptitle('Lawrence Stewart - Created using network().')
        degreeDistribution = nx.degree_histogram(G)
        print(degreeDistribution)
        plt.bar(np.arange(1,len(degreeDistribution)+1,1),[degreeDistribution[i]/len(G.nodes()) for i in range(len(degreeDistribution))],color='b')
        ax = plt.gca()
        ax.set_facecolor('#D9E6E8')
        plt.xticks(np.arange(1,len(degreeDistribution)+1,1))
        plt.title("Degree Distribution for Network Generated")
        plt.grid("on")
        plt.xlabel("Degree")
        plt.ylabel("Node Fraction")
        plt.show()

    return None

if __name__ == '__main__':
    #The code here should call analyze and generate the
    #figures that you are submitting with your code
    analyze(200,150)

reason: import numpy,import networkx


text:
from create_milp_dict import create_problem_dict
from uc_model import solve_milp, solution_to_schedule
from rl4uc.rl4uc.environment import make_env
import numpy as np
import pandas as pd

def test_kazarlis():
    """Kazarlis 1996 problem"""

    demand = np.array([700, 750, 850, 950, 1000, 1100,
                       1150, 1200, 1300, 1400, 1450, 1500, 1400,
                       1300, 1200, 1050, 1000, 1100, 1200,
                       1400, 1300, 1100, 900, 800])
    wind = np.array([0]*24)
    env_params = {'num_gen': 10, 'dispatch_freq_mins': 60}
    problem_dict = create_problem_dict(demand,
                                       wind,
                                       env_params,
                                       reserve_pct=10)
    solution = solve_milp(problem_dict)

    print("${:.2f}".format(solution.obj()))

#    assert (solution.obj() < 570000) and (solution.obj() > 550000)

    schedule = solution_to_schedule(solution, problem_dict)
    df = pd.DataFrame({'demand': demand, 'wind': wind})
    env = make_env(mode='test', profiles_df=df, **env_params)
    env.reset()
    cost = 0
    for a in schedule:
        print(a)
        o, r, d = env.step(a, deterministic=True)
        cost -= r
    print(cost)

if __name__ == '__main__':
    test_kazarlis()

max_stars_repo_name: pwdemars/pglib-uc

max_issues_repo_name: pwdemars/pglib-uc

max_forks_repo_name: pwdemars/pglib-uc


text:
import numpy as np
from matplotlib import pyplot as plt
import math
MAX_SPEED = 2
ACCELERATION = 0.5
DRAG = 0.3
TURN_SPEED=5

IMAGE = np.array([
[0,0,0,1,0,0,0],
[0,0,1,1,1,0,0],
[0,1,1,1,1,1,0],
[1,1,1,1,1,1,1],
[0,1,1,1,1,1,0],
[0,0,1,1,1,0,0],
[0,0,0,1,0,0,0]])

def main():
    position=(42 ,42)
    speed=0
    bearing=0
    acc=0
    turn=0

    plt.ion()
    fig, ax = plt.subplots()
    img = np.zeros((420,420))
    img[207:214,207:214]=IMAGE
    im = ax.imshow(img)

    for i in range(1000):
        acc+=np.random.rand()-0.5
        turn+=np.random.rand()-0.5
        acc=np.clip(acc,-1,1)
        turn=np.clip(turn,-1,1)
        (position,bearing,speed) = update_coords(position,bearing,speed,acc,turn)
        print(acc,turn)
        print(position)
        render(ax,im,position,bearing,speed)

def update_coords(position,bearing,speed,acceleration,turning):
    (x_pos,y_pos) = position
    speed = update_speed(speed,acceleration)
    bearing = (bearing + TURN_SPEED*turning) % 360
    x_pos += speed * math.sin(bearing*2*math.pi/360)
    y_pos += speed * math.cos(bearing*2*math.pi/360)
    return ((x_pos,y_pos),bearing,speed)

def update_speed(speed,acceleration):
    speed *= DRAG
    speed += acceleration*ACCELERATION
    speed = min(speed,MAX_SPEED) if speed > 0 else max(speed,-MAX_SPEED)
    return speed

def render(ax,im,position,bearing,speed):

    x_pos,y_pos = position
    img = np.zeros((420,420))
    x = int(x_pos*5)
    y = int(y_pos*5)
    img[x:x+7,y:y+7]=IMAGE

    # plt.scatter(x,y)
    # plt.show()
    im.set_data(img)
    ax.set_title(f"bearing : {bearing}, speed: {speed}")

    plt.pause(0.001)
    plt.draw()

if __name__ == "__main__":
    main()

max_stars_repo_name: AndrewPaulChester/sage-code

max_issues_repo_name: AndrewPaulChester/sage-code

max_forks_repo_name: AndrewPaulChester/sage-code