text:
# Copyright (c) 2012 The Chromium OS Authors. All rights reserved.
# Use of this source code is governed by a BSD-style license that can be
# found in the LICENSE file.
"""Interface to power measurement capabilities of TPS65090 PMIC.

Device contains 14 seperate ADCs allowing 12 power measurements of:
  1. Main AC input (via ADC VAC(0) & IAC(2))
  2. Battery Charging (via ADCs VBAT(1) & IBAT(3))
  3-5. Current measurement of 3 DC->DC converters.  Voltage must be provided.
  6-12. Current measurement of 7 power fets.  Voltage must be provided.

Module provides interface to register and read these power sources.

For example:

    i2c = TPS65090I2c(4)
    pwr = TPS65090Power(i2c)

    pwr_names = []
    pwr.register('ac', new_name='ac_mw')
    pwr_names.append('ac_mw')
    pwr.register('bat', new_name='bat_mw')
    pwr_names.append('bat_mw')

    voltages = [5, 3.3, 1.35, 40, 5, 3.3, 3.3, 3.3, 3.3, 5]
    for i,name in enumerate(AD_OUT_DEFAULT_NAMES[4:]):
        new_name = new_name='%s_mw' % name,
        pwr.register(name, new_name=name, voltage=voltages[i])
        pwr_names.append(new_name)

    for name in new_names:
        print '%s = %2.f\n' % (name, pwr.read_mw(name))

TODO(tbroch) Work with vendor to resolve discrepancies in measurements
"""
import commands
import logging
import numpy
import sys
import time

# TODO(tbroch) maybe this should be optarg
MAX_STATS = 20

# Lower 2bits configureable via Texas Instruments
I2C_SLV_DEFAULT = 0x48

REG_IDX_AD_CTRL = 0x16
# OUT1 is lower byte, OUT2 is upper byte of 10-bit ADC
REG_IDX_AD_OUT1 = 0x17
REG_IDX_AD_OUT2 = 0x18

# AD_CTRL POR value
# TODO(tbroch): Verify why <7> (reserved) != '1' despite datasheet POR
AD_CTRL_POR = 0x20

# Start a conversion measurement.  ~1mA of additional load to TPS65090 Vin
AD_CTRL_START = 0x40
# Latest conversion is completed
AD_CTRL_EOC   = 0x20
# ADC conversion vref is enabled.  ~.5mA of additional load to TPS65090 Vin
AD_CTRL_EN    = 0x10

AD_OUT_RANGE = 2**10 - 1
AD_OUT_LSB = [17./AD_OUT_RANGE,
              17./AD_OUT_RANGE,
              3./AD_OUT_RANGE,
              4./AD_OUT_RANGE,
              5./AD_OUT_RANGE,
              5./AD_OUT_RANGE,
              5./AD_OUT_RANGE,
              1.1/AD_OUT_RANGE,
              0.22/AD_OUT_RANGE,
              3.3/AD_OUT_RANGE,
              1.1/AD_OUT_RANGE,
              1.1/AD_OUT_RANGE,
              1.1/AD_OUT_RANGE,
              1.1/AD_OUT_RANGE]

AD_OUT_DEFAULT_NAMES = ['vac', 'vbat', 'ac', 'bat', 'dcdc1', 'dcdc2',
                        'dcdc3', 'fet1', 'fet2', 'fet3', 'fet4', 'fet5',
                        'fet6', 'fet7']

class TPS65090Error(Exception):
    """Exception class for TPS65090."""

class TPS65090I2c(object):
    """Class for TPS65090 I2c communication."""
    # -f to force, -y to respond 'yes' to interactive queries
    _GET_CMD = 'i2cget -f -y'
    _SET_CMD = 'i2cset -f -y'
    # TODO(tbroch) deprecate retry once there are no reports of errors on i2c
    _RETRY_MAX = 1

    def __init__(self, bus, slv=I2C_SLV_DEFAULT):
        """Constructor."""
        self._bus = bus
        self._slv = slv

    def _do_cmd(self, cmd):
        """Run command and return results.

        Args:
          cmd: string, i2c command to run

        Returns:
          string, output from shell command

        Raises:
          TPS65090Error: Bad status returned from shell command after retries
            exhausted

        """
        retry = self._RETRY_MAX
        status = 1
        while status and retry:
            (status, output) = commands.getstatusoutput(cmd)
            logging.debug('cmd: %s :: status: %d :: output %s',
                          cmd, status, output)
            retry -= 1
        if status:
            raise TPS65090Error('Non-zero status from i2c get returned '
                                'status:%d' % status)
        return output

    def get(self, index):
        """Get I2c value.

        Args:
          index: integer, register index to write byte to

        """
        cmd = '%s %d 0x%x 0x%x b' % (self._GET_CMD, self._bus, self._slv, index)
        output = self._do_cmd(cmd)
        return int(output, 0)

    def set(self, index, val):
        """Set I2c value.

        Args:
          index: integer, register index to write byte to
          val: integer, byte value to writ to index
        """
        cmd = '%s %d 0x%x 0x%x 0x%x' % \
            (self._SET_CMD, self._bus, self._slv, index, val)
        self._do_cmd(cmd)

class TPS65090ADC(object):
    """Class for TPS65090 ADC."""
    # Number of times to retry read of AD_CTRL for EOC
    _RETRY_MAX = 10

    def __init__(self, i2c, sensor_id):
        """."""
        self._i2c = i2c
        self._sensor_id = sensor_id

    def read_adc(self):
        """Read value for ADC.

        Steps required:
          1. Write ADC value into AD_CTRL[3:0]
          2. Enable measurement
          3. Poll conversion ready
          4. Read AD_OUT1
          5. Read AD_OUT2
          6. merge AD_OUT2[1:0] << 8 | AD_OUT1
          7. return 10bit * units_per_lsb

        Raises:
          TPS65090Error: if ADC conversion fails to complete
        """
        wr_val = AD_CTRL_EN | self._sensor_id
        self._i2c.set(REG_IDX_AD_CTRL, wr_val)
        wr_val |= AD_CTRL_START
        self._i2c.set(REG_IDX_AD_CTRL, wr_val)
        time.sleep(.01)

        rd_val = self._i2c.get(REG_IDX_AD_CTRL)
        retry = self._RETRY_MAX
        while retry and not (rd_val & AD_CTRL_EOC):
            time.sleep(.01)
            rd_val = self._i2c.get(REG_IDX_AD_CTRL)
            retry -= 1

        if not retry:
            raise TPS65090Error('Failed to see ADC conversion complete')

        rd_val1 = self._i2c.get(REG_IDX_AD_OUT1)
        rd_val2 = self._i2c.get(REG_IDX_AD_OUT2)
        adc_val = ((rd_val2 << 8) | rd_val1) * AD_OUT_LSB[self._sensor_id]
        logging.debug('adc%02d = %f', self._sensor_id, adc_val)
        return adc_val

class TPS65090Power(object):
    """Class for TPS65090 PMIC Power.

    Attributes:
      _i2c: object, TPS65090I2c instance to communicate to PMIC over I2c
      _voltage: list of voltages.  Either float or TPS65090ADC instances (for
        AC & battery)
      _current: list of currents.  All are TPS65090ADC instances
    """
    def __init__(self, i2c):
        """Constructor.

        Args:
          i2c: object, TPS65090I2c instance to communicate to PMIC over I2c
        """
        self._i2c = i2c
        self._voltage = {}
        self._current = {}

        val = self._i2c.get(REG_IDX_AD_CTRL)
        if val != AD_CTRL_POR:
            raise TPS65090Error('Unable to read POR value for AD_CTRL')
        self._i2c.set(REG_IDX_AD_CTRL, AD_CTRL_EN)

    def __del__(self):
        self._i2c.set(REG_IDX_AD_CTRL, AD_CTRL_POR)

    def register(self, default_current_name, new_name=None, voltage=None):
        """Register a power measurement.

        The TPS65090 provides two built-in power measurements by way
        of voltage and current measurements for the AC input and the
        battery charge.  The remaining 10 current sensors can be made
        to measure power with voltage provided.

        Args:
        default_current_name: string, current sensor default name per datasheet.
        new_name: string, name to call power measurement.  Defaults to
          default_current_name.
        voltage: float, value for potential of rail

        Raises:
          TPS65090Error: if voltage isn't provided for DCDC or FET sensors
        """
        sensor_id = AD_OUT_DEFAULT_NAMES.index(default_current_name)
        if not new_name:
            new_name = default_current_name

        if voltage is None:
            if sensor_id not in [2, 3]:
                raise TPS65090Error('Must provide nominal voltage')
            if sensor_id == 2:
                self._voltage[new_name] = TPS65090ADC(self._i2c, 0)
            else:
                self._voltage[new_name] = TPS65090ADC(self._i2c, 1)
        else:
            self._voltage[new_name] = float(voltage)

        self._current[new_name] = TPS65090ADC(self._i2c, sensor_id)
        logging.info('registered %s ( %s ) @ sensor %d', new_name,
                     default_current_name, sensor_id)

    def read_mw(self, name):
        """Read power in milliwatts.

        Returns:
          float, power in milliwatts

        Raises:
          TPS65090Error: if current ADC is unable to be found
        """
        if name not in self._current:
            raise TPS65090Error('Unable to measure current.  Did you '
                                'register the power?')

        if type(self._voltage[name]) is float:
            voltage = self._voltage[name]
        else:
            voltage = self._voltage[name].read_adc()
        current = self._current[name].read_adc()
        return current * voltage * 1000

def main():
    """main loop to read power."""
    logging.basicConfig(level=logging.INFO)
    logging.info('start')

    i2c = TPS65090I2c(4)
    pwr = TPS65090Power(i2c)

    pwr_names = []
    pwr.register('ac', new_name='ac_mw')
    pwr_names.append('ac_mw')
    pwr.register('bat', new_name='bat_mw')
    pwr_names.append('bat_mw')

    voltages = [5, 3.3, 1.35, 12, 5, 3.3, 3.3, 3.3, 3.3, 5]
    names = ['pp5000_aux', 'pp3300_aux', 'pp1350_ddr', 'bl', 'pp5000_video',
              'pp3300_wwan', 'pp3300_sdcard', 'pp3300_cam', 'pp3300_lcd',
              'pp5000_ts']

    for i, name in enumerate(AD_OUT_DEFAULT_NAMES[4:]):
        new_name = '%s_mw' % names[i]
        pwr.register(name, new_name=new_name, voltage=voltages[i])
        pwr_names.append(new_name)

    rails = pwr_names
    if len(sys.argv) > 1:
        rails = sys.argv[1:]

    stats = {}
    for name in rails:
        stats[name] = numpy.array([])

    while True:
        for name in rails:
            mw_val = pwr.read_mw(name)
            if len(stats[name]) == MAX_STATS:
                stats[name] = numpy.delete(stats[name], 0)
            stats[name] = numpy.append(stats[name], mw_val)
            print '%s:%.f' % (name, stats[name].mean())
        # delay in order to not raise power via the script
        # TODO(tbroch) sleep should be configureable
        time.sleep(2)

if __name__ == '__main__':
    try:
        main()
    except KeyboardInterrupt:
        sys.exit(0)

max_stars_repo_name: mmind/servo-hdctools

max_issues_repo_name: mmind/servo-hdctools

max_forks_repo_name: mmind/servo-hdctools


text:
#!/usr/bin/env python3
"""Script to train the deterministic supervised adversarial autoencoder."""
from pathlib import Path
import random as rn
import time

import joblib
from sklearn.preprocessing import RobustScaler, OneHotEncoder
import numpy as np
import tensorflow as tf

from utils import COLUMNS_NAME, load_dataset
from models import make_encoder_model_v1, make_decoder_model_v1, make_discriminator_model_v1

PROJECT_ROOT = Path.cwd()

def main():
    """Train the normative method on the bootstrapped samples.

    The script also the scaler and the demographic data encoder.
    """
    # ----------------------------------------------------------------------------
    n_bootstrap = 1000
    model_name = 'supervised_aae'

    participants_path = PROJECT_ROOT / 'data' / 'BIOBANK' / 'participants.tsv'
    freesurfer_path = PROJECT_ROOT / 'data' / 'BIOBANK' / 'freesurferData.csv'
    # ----------------------------------------------------------------------------
    bootstrap_dir = PROJECT_ROOT / 'outputs' / 'bootstrap_analysis'
    ids_dir = bootstrap_dir / 'ids'

    model_dir = bootstrap_dir / model_name
    model_dir.mkdir(exist_ok=True)

    # ----------------------------------------------------------------------------
    # Set random seed
    random_seed = 42
    tf.random.set_seed(random_seed)
    np.random.seed(random_seed)
    rn.seed(random_seed)

    for i_bootstrap in range(n_bootstrap):
        ids_filename = 'cleaned_bootstrap_{:03d}.csv'.format(i_bootstrap)
        ids_path = ids_dir / ids_filename

        bootstrap_model_dir = model_dir / '{:03d}'.format(i_bootstrap)
        bootstrap_model_dir.mkdir(exist_ok=True)

        # ----------------------------------------------------------------------------
        # Loading data
        dataset_df = load_dataset(participants_path, ids_path, freesurfer_path)

        # ----------------------------------------------------------------------------
        x_data = dataset_df[COLUMNS_NAME].values

        tiv = dataset_df['EstimatedTotalIntraCranialVol'].values
        tiv = tiv[:, np.newaxis]

        x_data = (np.true_divide(x_data, tiv)).astype('float32')

        scaler = RobustScaler()
        x_data_normalized = scaler.fit_transform(x_data)

        # ----------------------------------------------------------------------------
        age = dataset_df['Age'].values[:, np.newaxis].astype('float32')
        enc_age = OneHotEncoder(sparse=False)
        one_hot_age = enc_age.fit_transform(age)

        gender = dataset_df['Gender'].values[:, np.newaxis].astype('float32')
        enc_gender = OneHotEncoder(sparse=False)
        one_hot_gender = enc_gender.fit_transform(gender)

        y_data = np.concatenate((one_hot_age, one_hot_gender), axis=1).astype('float32')

        # -------------------------------------------------------------------------------------------------------------
        # Create the dataset iterator
        batch_size = 256
        n_samples = x_data.shape[0]

        train_dataset = tf.data.Dataset.from_tensor_slices((x_data_normalized, y_data))
        train_dataset = train_dataset.shuffle(buffer_size=n_samples)
        train_dataset = train_dataset.batch(batch_size)

        # -------------------------------------------------------------------------------------------------------------
        # Create models
        n_features = x_data_normalized.shape[1]
        n_labels = y_data.shape[1]
        h_dim = [100, 100]
        z_dim = 20

        encoder = make_encoder_model_v1(n_features, h_dim, z_dim)
        decoder = make_decoder_model_v1(z_dim + n_labels, n_features, h_dim)
        discriminator = make_discriminator_model_v1(z_dim, h_dim)

        # -------------------------------------------------------------------------------------------------------------
        # Define loss functions
        cross_entropy = tf.keras.losses.BinaryCrossentropy(from_logits=True)
        mse = tf.keras.losses.MeanSquaredError()
        accuracy = tf.keras.metrics.BinaryAccuracy()

        def discriminator_loss(real_output, fake_output):
            loss_real = cross_entropy(tf.ones_like(real_output), real_output)
            loss_fake = cross_entropy(tf.zeros_like(fake_output), fake_output)
            return loss_fake + loss_real

        def generator_loss(fake_output):
            return cross_entropy(tf.ones_like(fake_output), fake_output)

        # -------------------------------------------------------------------------------------------------------------
        # Define optimizers
        base_lr = 0.0001
        max_lr = 0.005

        step_size = 2 * np.ceil(n_samples / batch_size)

        ae_optimizer = tf.keras.optimizers.Adam(lr=base_lr)
        dc_optimizer = tf.keras.optimizers.Adam(lr=base_lr)
        gen_optimizer = tf.keras.optimizers.Adam(lr=base_lr)

        # -------------------------------------------------------------------------------------------------------------
        # Training function
        @tf.function
        def train_step(batch_x, batch_y):
            # -------------------------------------------------------------------------------------------------------------
            # Autoencoder
            with tf.GradientTape() as ae_tape:
                encoder_output = encoder(batch_x, training=True)
                decoder_output = decoder(tf.concat([encoder_output, batch_y], axis=1), training=True)

                # Autoencoder loss
                ae_loss = mse(batch_x, decoder_output)

            ae_grads = ae_tape.gradient(ae_loss, encoder.trainable_variables + decoder.trainable_variables)
            ae_optimizer.apply_gradients(zip(ae_grads, encoder.trainable_variables + decoder.trainable_variables))

            # -------------------------------------------------------------------------------------------------------------
            # Discriminator
            with tf.GradientTape() as dc_tape:
                real_distribution = tf.random.normal([batch_x.shape[0], z_dim], mean=0.0, stddev=1.0)
                encoder_output = encoder(batch_x, training=True)

                dc_real = discriminator(real_distribution, training=True)
                dc_fake = discriminator(encoder_output, training=True)

                # Discriminator Loss
                dc_loss = discriminator_loss(dc_real, dc_fake)

                # Discriminator Acc
                dc_acc = accuracy(tf.concat([tf.ones_like(dc_real), tf.zeros_like(dc_fake)], axis=0),
                                  tf.concat([dc_real, dc_fake], axis=0))

            dc_grads = dc_tape.gradient(dc_loss, discriminator.trainable_variables)
            dc_optimizer.apply_gradients(zip(dc_grads, discriminator.trainable_variables))

            # -------------------------------------------------------------------------------------------------------------
            # Generator (Encoder)
            with tf.GradientTape() as gen_tape:
                encoder_output = encoder(batch_x, training=True)
                dc_fake = discriminator(encoder_output, training=True)

                # Generator loss
                gen_loss = generator_loss(dc_fake)

            gen_grads = gen_tape.gradient(gen_loss, encoder.trainable_variables)
            gen_optimizer.apply_gradients(zip(gen_grads, encoder.trainable_variables))

            return ae_loss, dc_loss, dc_acc, gen_loss

        # -------------------------------------------------------------------------------------------------------------
        # Training loop
        global_step = 0
        n_epochs = 200
        gamma = 0.98
        scale_fn = lambda x: gamma ** x
        for epoch in range(n_epochs):
            start = time.time()

            epoch_ae_loss_avg = tf.metrics.Mean()
            epoch_dc_loss_avg = tf.metrics.Mean()
            epoch_dc_acc_avg = tf.metrics.Mean()
            epoch_gen_loss_avg = tf.metrics.Mean()

            for _, (batch_x, batch_y) in enumerate(train_dataset):
                global_step = global_step + 1
                cycle = np.floor(1 + global_step / (2 * step_size))
                x_lr = np.abs(global_step / step_size - 2 * cycle + 1)
                clr = base_lr + (max_lr - base_lr) * max(0, 1 - x_lr) * scale_fn(cycle)
                ae_optimizer.lr = clr
                dc_optimizer.lr = clr
                gen_optimizer.lr = clr

                ae_loss, dc_loss, dc_acc, gen_loss = train_step(batch_x, batch_y)

                epoch_ae_loss_avg(ae_loss)
                epoch_dc_loss_avg(dc_loss)
                epoch_dc_acc_avg(dc_acc)
                epoch_gen_loss_avg(gen_loss)

            epoch_time = time.time() - start

            print('{:4d}: TIME: {:.2f} ETA: {:.2f} AE_LOSS: {:.4f} DC_LOSS: {:.4f} DC_ACC: {:.4f} GEN_LOSS: {:.4f}' \
                  .format(epoch, epoch_time,
                          epoch_time * (n_epochs - epoch),
                          epoch_ae_loss_avg.result(),
                          epoch_dc_loss_avg.result(),
                          epoch_dc_acc_avg.result(),
                          epoch_gen_loss_avg.result()))

        # Save models
        encoder.save(bootstrap_model_dir / 'encoder.h5')
        decoder.save(bootstrap_model_dir / 'decoder.h5')
        discriminator.save(bootstrap_model_dir / 'discriminator.h5')

        # Save scaler
        joblib.dump(scaler, bootstrap_model_dir / 'scaler.joblib')

        joblib.dump(enc_age, bootstrap_model_dir / 'age_encoder.joblib')
        joblib.dump(enc_gender, bootstrap_model_dir / 'gender_encoder.joblib')

if __name__ == "__main__":
    main()

ext: py

reason: import numpy


text:
"""Affine connections."""

import autograd
from scipy.optimize import minimize

import geomstats.backend as gs
from geomstats.integrator import integrate

N_STEPS = 10
EPSILON = 1e-3

class Connection(object):
    """Classe for affine connections."""

    def __init__(self, dimension):
        self.dimension = dimension

    def christoffels(self, base_point):
        """Christoffel symbols associated with the connection.

        Parameters
        ----------
        base_point : array-like, shape=[n_samples, dimension]

        Returns
        -------
        gamma : array-like, shape=[n_samples, dimension, dimension, dimension]
        """
        raise NotImplementedError(
            'The Christoffel symbols are not implemented.')

    def connection(self, tangent_vector_a, tangent_vector_b, base_point):
        """Covariant derivative.

        Connection applied to `tangent_vector_b` in the direction of
        `tangent_vector_a`, both tangent at `base_point`.

        Parameters
        ----------
        tangent_vec_a : array-like, shape=[n_samples, dimension]
        tangent_vec_b : array-like, shape=[n_samples, dimension]
        base_point : array-like, shape=[n_samples, dimension]
        """
        raise NotImplementedError(
            'connection is not implemented.')

    def geodesic_equation(self, position, velocity):
        """Compute the geodesic ODE associated with the connection.

        Parameters
        ----------
        velocity : array-like, shape=[n_samples, dimension]
        position : array-like, shape=[n_samples, dimension]
            the position at which to compute the geodesic ODE

        Returns
        -------
        geodesic_ode : array-like, shape=[n_samples, dimension]
            value of the vector field to be integrated at position
        """
        gamma = self.christoffels(position)
        return - gs.einsum('...kij,...i,...j->...k', gamma, velocity,
                           velocity)

    def exp(self, tangent_vec, base_point, n_steps=N_STEPS, step='euler'):
        """Exponential map associated to the affine connection.

        Exponential map at base_point of tangent_vec computed by integration
        of the geodesic equation (initial value problem), using the
        christoffel symbols

        Parameters
        ----------
        tangent_vec : array-like, shape=[n_samples, dimension]
        base_point : array-like, shape=[n_samples, dimension]
        n_steps : int
            the number of discrete time steps to take in the integration
        step : str, {'euler', 'rk4'}
            the numerical scheme to use for integration

        Returns
        -------
        exp : array-like, shape=[n_samples, dimension]
        """
        tangent_vec = gs.to_ndarray(tangent_vec, to_ndim=2)
        base_point = gs.to_ndarray(base_point, to_ndim=2)
        initial_state = (base_point, tangent_vec)
        flow, _ = integrate(
            self.geodesic_equation, initial_state, n_steps=n_steps, step=step)
        return flow[-1]

    def log(self, point, base_point, n_steps=N_STEPS, step='euler'):
        """Compute logarithm map associated to the affine connection.

        Solve the boundary value problem associated to the geodesic equation
        using the Christoffel symbols and conjugate gradient descent.

        Parameters
        ----------
        point : array-like, shape=[n_samples, dimension]
        base_point : array-like, shape=[n_samples, dimension]
        n_steps : int
            the number of discrete time steps to take in the integration
        step : str, {'euler', 'rk4'}
            the numerical scheme to use for integration

        Returns
        -------
        tangent_vec : array-like, shape=[n_samples, dimension]
        """
        point = gs.to_ndarray(point, to_ndim=2)
        base_point = gs.to_ndarray(base_point, to_ndim=2)
        n_samples = len(base_point)

        def objective(velocity):
            """Define the objective function."""
            velocity = velocity.reshape(base_point.shape)
            delta = self.exp(velocity, base_point, n_steps, step) - point
            loss = 1. / self.dimension * gs.sum(delta ** 2, axis=1)
            return 1. / n_samples * gs.sum(loss)

        objective_grad = autograd.elementwise_grad(objective)

        def objective_with_grad(velocity):
            """Create helpful objective func wrapper for autograd comp."""
            return objective(velocity), objective_grad(velocity)

        tangent_vec = gs.random.rand(base_point.size)
        res = minimize(
            objective_with_grad, tangent_vec, method='L-BFGS-B', jac=True,
            options={'disp': False, 'maxiter': 25})

        tangent_vec = res.x
        tangent_vec = gs.reshape(tangent_vec, base_point.shape)
        return tangent_vec

    def pole_ladder_step(self, base_point, next_point, base_shoot):
        """Compute one Pole Ladder step.

        One step of pole ladder scheme [LP2013a]_ using the geodesic to
        transport along as diagonal of the parallelogram.

        Parameters
        ----------
        base_point : array-like
            shape=[n_samples, dimension] or shape=[1, dimension]
        next_point : array-like
            shape=[n_samples, dimension] or shape=[1, dimension]
        base_shoot : array-like
            shape=[n_samples, dimension] or shape=[1, dimension]

        Returns
        -------
        transported_tangent_vector : array-like
            shape=[n_samples, dimension] or shape=[1, dimension]
        end_point : array-like
            shape=[n_samples, dimension] or shape=[1, dimension]

        References
        ----------
        .. [LP2013a] Marco Lorenzi, Xavier Pennec. Efficient Parallel Transport
         of Deformations in Time Series of Images: from Schild's to
         Pole Ladder. Journal of Mathematical Imaging and Vision, Springer
         Verlag, 2013,50 (1-2), pp.5-17. ⟨10.1007/s10851-013-0470-3⟩
        """
        mid_tangent_vector_to_shoot = 1. / 2. * self.log(
            base_point=base_point,
            point=next_point)

        mid_point = self.exp(
            base_point=base_point,
            tangent_vec=mid_tangent_vector_to_shoot)

        tangent_vector_to_shoot = - self.log(
            base_point=mid_point,
            point=base_shoot)

        end_shoot = self.exp(
            base_point=mid_point,
            tangent_vec=tangent_vector_to_shoot)

        transported_tangent_vector = - self.log(
            base_point=next_point, point=end_shoot)

        end_point = self.exp(
            base_point=next_point,
            tangent_vec=transported_tangent_vector)

        return transported_tangent_vector, end_point

    def pole_ladder_parallel_transport(
            self, tangent_vec_a, tangent_vec_b, base_point, n_steps=1):
        """Approximate parallel transport using the pole ladder scheme.

        Approximate Parallel transport using the pole ladder scheme [LP2013b]_
        [GJSP2019]_. `tangent_vec_a` is transported along the geodesic starting
        at the base_point with initial tangent vector `tangent_vec_b`.

        Returns a tangent vector at the point
        exp_(`base_point`)(`tangent_vec_b`).

        Parameters
        ----------
        tangent_vec_a : array-like, shape=[n_samples, dimension]
                                   or shape=[1, dimension]
        tangent_vec_b : array-like, shape=[n_samples, dimension]
                                   or shape=[1, dimension]
        base_point : array-like, shape=[n_samples, dimension]
                                or shape=[1, dimension]
        n_steps: int
            the number of pole ladder steps

        Returns
        -------
        transported_tangent_vector : array-like, shape=[n_samples, dimension]
                                                or shape=[1, dimension]

        References
        ----------
        .. [LP2013b] Marco Lorenzi, Xavier Pennec. Efficient Parallel Transpor
          of Deformations in Time Series of Images: from Schild's to
          Pole Ladder.Journal of Mathematical Imaging and Vision, Springer
          Verlag, 2013, 50 (1-2), pp.5-17. ⟨10.1007/s10851-013-0470-3⟩

        .. [GJSP2019] N. Guigui, Shuman Jia, Maxime Sermesant, Xavier Pennec.
          Symmetric Algorithmic Components for Shape Analysis with
          Diffeomorphisms. GSI 2019, Aug 2019, Toulouse, France. pp.10.
          ⟨hal-02148832⟩
        """
        current_point = gs.copy(base_point)
        transported_tangent_vector = gs.copy(tangent_vec_a)
        base_shoot = self.exp(base_point=current_point,
                              tangent_vec=transported_tangent_vector)
        for i_point in range(0, n_steps):
            frac_tangent_vector_b = (i_point + 1) / n_steps * tangent_vec_b
            next_point = self.exp(
                base_point=base_point,
                tangent_vec=frac_tangent_vector_b)
            transported_tangent_vector, base_shoot = self.pole_ladder_step(
                base_point=current_point,
                next_point=next_point,
                base_shoot=base_shoot)
            current_point = next_point

        return transported_tangent_vector

    def riemannian_curvature(self, base_point):
        """Compute Riemannian curvature tensor associated with the connection.

        Parameters
        ----------
        base_point: array-like, shape=[n_samples, dimension]
                                or shape=[1, dimension]
        """
        raise NotImplementedError(
            'The Riemannian curvature tensor is not implemented.')

    def geodesic(self, initial_point,
                 end_point=None, initial_tangent_vec=None,
                 point_type='vector'):
        """Generate parameterized function for the geodesic curve.

        Geodesic curve defined by either:
        - an initial point and an initial tangent vector,
        - an initial point and an end point.

        Parameters
        ----------
        initial_point
        end_point
        initial_tangent_vec
        point_type

        Returns
        -------
        path : callable
            the time parameterized geodesic curve.
        """
        point_ndim = 1
        if point_type == 'matrix':
            point_ndim = 2

        initial_point = gs.to_ndarray(
            initial_point, to_ndim=point_ndim + 1)

        if end_point is None and initial_tangent_vec is None:
            raise ValueError('Specify an end point or an initial tangent '
                             'vector to define the geodesic.')
        if end_point is not None:
            end_point = gs.to_ndarray(
                end_point, to_ndim=point_ndim + 1)
            shooting_tangent_vec = self.log(point=end_point,
                                            base_point=initial_point)
            if initial_tangent_vec is not None:
                assert gs.allclose(shooting_tangent_vec, initial_tangent_vec)
            initial_tangent_vec = shooting_tangent_vec
        initial_tangent_vec = gs.array(initial_tangent_vec)
        initial_tangent_vec = gs.to_ndarray(
            initial_tangent_vec, to_ndim=point_ndim + 1)

        def path(t):
            """Generate parameterized function for geodesic curve."""
            t = gs.cast(t, gs.float32)
            t = gs.to_ndarray(t, to_ndim=1)
            t = gs.to_ndarray(t, to_ndim=2, axis=1)
            new_initial_point = gs.to_ndarray(
                initial_point,
                to_ndim=point_ndim + 1)
            new_initial_tangent_vec = gs.to_ndarray(
                initial_tangent_vec,
                to_ndim=point_ndim + 1)

            if point_type == 'vector':
                tangent_vecs = gs.einsum(
                    'il,nk->ik',
                    t,
                    new_initial_tangent_vec)
            elif point_type == 'matrix':
                tangent_vecs = gs.einsum(
                    'il,nkm->ikm',
                    t,
                    new_initial_tangent_vec)

            point_at_time_t = self.exp(
                tangent_vec=tangent_vecs,
                base_point=new_initial_point)
            return point_at_time_t

        return path

    def torsion(self, base_point):
        """Compute torsion tensor associated with the connection.

        Parameters
        ----------
        base_point: array-like, shape=[n_samples, dimension]
                                or shape=[1, dimension]
        """
        raise NotImplementedError(
            'The torsion tensor is not implemented.')

reason: from scipy


text:
import gym
from gym import error, spaces, utils
from gym.utils import seeding
import math
import numpy as np
from gym.envs.classic_control import rendering

class NavigateEnv(gym.Env):
  metadata = {
  	'render.modes': ['human', 'rgb_array'],
  	'video.frames_per_second' : 28
  }

  def __init__(self):

    # Bot dynamics
    self.min_position = np.array([0, 0])
    self.max_position = np.array([1630 / 5, 1948 / 5])
    self.goal_position = np.array([0, 1066 / 5])
    self.start_position = np.array([1630 / 5, 1066 / 5])
    self.min_step = 0
    self.max_step = 2.0
    self.frames_per_step = 4
    self.distance_to_goal = 0.0
    self.prev_distance_to_goal = 0.0

    # If below -1 reward for n steps in a row, terminate game
    self.fault = [0,5] # ---> [current count, n]

    # Viewer
    self.screen_width = 700
    self.screen_height = 700
    self.grid_boxes = 100
    self.bot_width = 40
    self.bot_height = 20

    # Spaces
    self.action_space = spaces.Box(
  		low=np.array([-self.max_step, -self.max_step]),
  		high=np.array([self.max_step, self.max_step]))
    self.observation_space = spaces.Box(
    	low=np.array([self.min_position[0], self.min_position[1]]),
    	high=np.array([self.max_position[0], self.max_position[1]]))

    self.viewer = None
    self.state = None

    self.reset()

  def _convert_graphics_to_viewer_coords(graphics_coords):
    return (graphics_coords[1], graphics_coords[0])

  def step(self, action):
    assert self.action_space.contains(action), "%r (%s) invalid" % (action, type(action))

    position = self.state
    position = (position[0] + action[0], position[1] + action[1])

    if (position[0] > self.max_position[0]): position = (self.max_position[0], position[1])
    if (position[1] > self.max_position[1]): position = (position[0], self.max_position[1])
    if (position[0] < self.min_position[0]): position = (self.min_position[0], position[1])
    if (position[1] < self.min_position[1]): position = (position[0, self.min_position[1]])

    done = bool((position[0] == self.goal_position[0]
      and position[1] == self.goal_position[1])
      or self.fault[0] >= self.fault[1])

    # calculate reward
    reward = 0
    if (not done):

      # future considerations
      #  - penalize each step
      #  - give high final reward upon reaching goal
      #  - increase reward inverse proportionally to gap between the agent and its goal

      self.distance_to_goal = np.linalg.norm(position - self.goal_position)
      if (self.distance_to_goal > self.prev_distance_to_goal):
        reward = -1.0
        self.fault[0] += 1
      elif (self.distance_to_goal == self.prev_distance_to_goal):
        reward = 0.0
      else:
        reward = 1.0
        self.fault[0] = 0

    self.state = position
    return self.state, reward, done, {}

  def reset(self):
    self.state = self.start_position
    return np.array(self.state)

  def render(self, mode='human'):

    if self.viewer is None:
      self.viewer = rendering.Viewer(self.screen_width, self.screen_height)

      # # list of evenly spaced x and y coordinates for the grid boxes
      # xs = np.linspace(self.min_position[0], self.max_position[1], self.grid_boxes)
      # ys = self._height(xs)
      # xys = list(zip((xs-self.min_position[0])*self.scale, ys*self.scale))

      # nikhil - may use later when creating a line that trails the bot to show the entire path taken
      # track = rendering.make_polyline(xys)
      # track.set_linewidth(4)
      # self.viewer.add_geom(track)

      # bot
      l = self.start_position[0] - self.bot_width / 2
      r = self.start_position[0] + self.bot_width / 2
      t = self.start_position[1] + self.bot_width / 2
      b = self.start_position[1] - self.bot_width / 2
      bot = rendering.FilledPolygon([(l,b), (l,t), (r,t), (r,b)])
      bot.add_attr(rendering.Transform())
      self.bot_trans = rendering.Transform()
      bot.add_attr(self.bot_trans)
      self.viewer.add_geom(bot)

      # GOOOOOOAAAAAAAALLLLLLL #TODO
      flagx = self.goal_position[0] + 40
      flagy1 = self.goal_position[1] - 20
      flagy2 = self.goal_position[1] + 20

      flagpole = rendering.Line((flagx, flagy1), (flagx, flagy2))
      self.viewer.add_geom(flagpole)
      #flag = rendering.FilledPolygon([(flagx, flagy2), (flagx, flagy2-10), (flagx+25, flagy2-5)])
      #flag.set_color(.8,.8,0)
      #self.viewer.add_geom(flag)

      pos = self.state[0]

      # nikhil - may use later to show rotational movement of the bot
      #self.bot_trans.set_translation((np.linalg.norm(pos-self.min_position[0]))*self.scale, self._height(pos[0])*self.scale) #TODO
      #self.bot_trans.set_rotation(math.cos(3 * pos[0])) #TODO

      return self.viewer.render(return_rgb_array = mode == 'rgb_array')

  def close(self):
    if self.viewer: self.viewer.close()


text:
import numpy as np
import matplotlib.pyplot as plt
from sim21.data.chemsep import pure
from sim21.data.chemsep_consts import GAS_CONSTANT
from sim21.provider.cubic import PengRobinson
from sim21.provider.flash.basic import basic_flash_temp_press_2phase
from sim21.provider.flash.io import flash_press_vap_frac_2phase, flash_temp_vap_frac_2phase, flash_press_prop_2phase, \
    flash_temp_prop_2phase, flash_prop_vap_frac_2phase

def press_vs_rho():
    components = pure(['ETHANE', 'N-HEPTANE'])
    prov = PengRobinson(components)
    n = np.array([0.5, 0.5])

    press_range = np.arange(0.01, 10, 0.01)
    temp_range = [420, 500]

    for temp in temp_range:
        rho_vapor, rho_liquid = [], []
        norm_hres_vapor, norm_hres_liquid = [], []
        log_phi_vapor_C2, log_phi_liquid_C2 = [], []
        log_phi_vapor_C7, log_phi_liquid_C7 = [], []

        for press in press_range:
            ph_vap = prov.phase(temp, press*1e6, n, 'vap')
            ph_liq = prov.phase(temp, press*1e6, n, 'liq')

            rho_vapor.append(ph_vap.dens_mole)
            rho_liquid.append(ph_liq.dens_mole)

            norm_hres_vapor.append(-ph_vap.res_enthalpy_mole/(GAS_CONSTANT*temp))
            norm_hres_liquid.append(-ph_liq.res_enthalpy_mole/(GAS_CONSTANT*temp))

            log_phi_vapor_C2.append(ph_vap.log_phi[0])
            log_phi_vapor_C7.append(ph_vap.log_phi[1])

            log_phi_liquid_C2.append(ph_liq.log_phi[0])
            log_phi_liquid_C7.append(ph_liq.log_phi[1])

        plt.plot(rho_vapor, press_range)
        plt.plot(rho_liquid, press_range)
        plt.legend(['rho_vapor', 'rho_liquid'], loc='upper right')
        plt.title('Equimolar ethane-n_heptane mixture @ '+str(temp)+' K')
        plt.xlabel('Density (kg/m3)')
        plt.ylabel('Press (MPa)')
        plt.show()
        plt.close()

        plt.plot(press_range, norm_hres_vapor)
        plt.plot(press_range, norm_hres_liquid)
        plt.legend(['norm_hres_vapor', 'norm_hres_liquid'], loc='upper right')
        plt.title('Equimolar ethane-n_heptane mixture @ '+str(temp)+' K')
        plt.xlabel('Press (MPa)')
        plt.ylabel('-DH_res/RT')
        plt.show()
        plt.close()

        plt.plot(press_range, log_phi_vapor_C2)
        plt.plot(press_range, log_phi_vapor_C7)
        plt.plot(press_range, log_phi_liquid_C2)
        plt.plot(press_range, log_phi_liquid_C7)
        plt.legend(['log_phi_vapor_C2', 'log_phi_vapor_C7', 'log_phi_liquid_C2', 'log_phi_liquid_C7'], loc='upper right')
        plt.title('Equimolar ethane-n_heptane mixture @ ' + str(temp) + ' K')
        plt.xlabel('Press (MPa)')
        plt.ylabel('log_phi')
        plt.show()
        plt.close()

def rhos_vs_temp():
    components = pure(['ETHANE', 'N-HEPTANE'])
    prov = PengRobinson(components)
    comp = np.array([0.5, 0.5])
    temp_list = np.arange(300, 500)
    rho_vapor, rho_liquid = [], []
    press = 2.0*1e6

    for temp in temp_list:
        ph_vap = prov.phase(temp, press, comp, 'vap')
        ph_liq = prov.phase(temp, press, comp, 'liq')
        rho_vapor.append(ph_vap.dens_mole)
        rho_liquid.append(ph_liq.dens_mole)

    plt.plot(temp_list, rho_vapor)
    plt.plot(temp_list, rho_liquid)
    plt.legend(['rho_vapor', 'rho_liquid'], loc='upper right')
    plt.title('Equimolar ethane-n_heptane mixture')
    plt.xlabel('Temp. (K)')
    plt.ylabel('Density (kg/m3)')
    plt.show()
    plt.close()

def test_flash():
    components = pure(['ETHANE', 'N-BUTANE', 'N-HEPTANE'])
    prov = PengRobinson(components)
    temp, press = 300, 101325 * 5
    feed_comp = np.array([0.2, 0.3, 0.5])

    # Repeat flash
    temp = 406.315
    results = basic_flash_temp_press_2phase(prov, temp, press, feed_comp, prov.all_valid_components)
    print('temp:', temp, 'V/F:', results.phases, results.vap_frac_mole, results.vap.enthalpy_mole, results.enthalpy_mole)

    temp = 280.61
    results = basic_flash_temp_press_2phase(prov, temp, press, feed_comp, prov.all_valid_components)
    print('temp:', temp, 'V/F:', results.phases, results.vap_frac_mole)

    return results

def test_cavett_feed():
    components = pure(['NITROGEN', 'CARBON DIOXIDE', 'HYDROGEN SULFIDE',
                       'METHANE', 'ETHANE', 'PROPANE', 'N-BUTANE', 'ISOBUTANE',
                       'N-PENTANE', 'ISOPENTANE', 'N-HEXANE', 'N-HEPTANE',
                       'N-OCTANE', 'N-NONANE', 'N-DECANE',
                       'N-UNDECANE'])

    prov = PengRobinson(components)
    temp, press = 322.04, 0.8*1e6
    feed_comp = np.array([0.0131, 0.1816, 0.0124, 0.1096, 0.0876,
                          0.0838, 0.0563, 0.0221, 0.0413, 0.0289,
                          0.0645, 0.0953, 0.0675, 0.0610, 0.0304,
                          0.0444])

    result = None
    for temp in range(150, 500):
        results = basic_flash_temp_press_2phase(prov, temp, press, feed_comp, prov.all_valid_components)
        print('temp:', temp, 'V/F:', results.phases, results.vap_frac_mole)

    return results

def test_flash_press_vap_frac():
    components = pure(['ETHANE', 'N-BUTANE', 'N-HEPTANE'])
    prov = PengRobinson(components)
    temp, press = 300, 101325 * 12
    feed_comp = np.array([0.2, 0.3, 0.5])
    for vf in np.arange(0, 1.01, 0.01):
        r = flash_press_vap_frac_2phase(prov, press, vf, feed_comp)
        print('temp:', r.temp, 'press:', r.press, 'vap_frac', r.vap_frac_mole)

    components = pure(['NITROGEN', 'CARBON DIOXIDE', 'HYDROGEN SULFIDE',
                       'METHANE', 'ETHANE', 'PROPANE', 'N-BUTANE', 'ISOBUTANE',
                       'N-PENTANE', 'ISOPENTANE', 'N-HEXANE', 'N-HEPTANE',
                       'N-OCTANE', 'N-NONANE', 'N-DECANE',
                       'N-UNDECANE'])

    prov = PengRobinson(components)
    temp, press = 322.04, 0.5*1e6
    feed_comp = np.array([0.0131, 0.1816, 0.0124, 0.1096, 0.0876,
                          0.0838, 0.0563, 0.0221, 0.0413, 0.0289,
                          0.0645, 0.0953, 0.0675, 0.0610, 0.0304,
                          0.0444])
    for vf in np.arange(0, 1.01, 0.01):
        r = flash_press_vap_frac_2phase(prov, press, vf, feed_comp)
        print('temp:', r.temp, 'press:', r.press, 'vap_frac', r.vap_frac_mole)

def test_flash_temp_vap_frac():
    components = pure(['NITROGEN', 'CARBON DIOXIDE', 'HYDROGEN SULFIDE',
                       'METHANE', 'ETHANE', 'PROPANE', 'N-BUTANE', 'ISOBUTANE',
                       'N-PENTANE', 'ISOPENTANE', 'N-HEXANE', 'N-HEPTANE',
                       'N-OCTANE', 'N-NONANE', 'N-DECANE',
                       'N-UNDECANE'])

    prov = PengRobinson(components)
    temp, press = 316.5395439844929, 0.5*1e6
    feed_comp = np.array([0.0131, 0.1816, 0.0124, 0.1096, 0.0876,
                          0.0838, 0.0563, 0.0221, 0.0413, 0.0289,
                          0.0645, 0.0953, 0.0675, 0.0610, 0.0304,
                          0.0444])

    vf = 0.5
    # This flash can be a little more buggy ... we will fix later
    r = flash_temp_vap_frac_2phase(prov, temp, vf, feed_comp)
    print('temp:', r.temp, 'press:', r.press, 'vap_frac', r.vap_frac_mole)

def test_flash_press_prop():
    components = pure(['ETHANE', 'N-BUTANE', 'N-HEPTANE'])
    prov = PengRobinson(components)
    temp, press = 300, 101325 * 5
    feed_comp = np.array([0.2, 0.3, 0.5])

    # Repeat flash
    for temp in [50, 100, 150, 200, 300, 320, 330, 380, 400, 480, 500, 600, 800, 1000]:
        results = basic_flash_temp_press_2phase(prov, temp, press, feed_comp, prov.all_valid_components)
        print('temp:', temp, 'V/F:', results.vap_frac_mole)

        results = flash_press_prop_2phase(prov, press, 'entropy_mole', results.entropy_mole, 0, feed_comp, valid=None, previous=None)
        print('temp:', results.temp, 'V/F:', results.vap_frac_mole)
        print('---')

    # Repeat flash
    for temp in [50, 100, 150, 200, 300, 320, 330, 380, 400, 480, 500, 600, 800, 1000]:
        results = basic_flash_temp_press_2phase(prov, temp, press, feed_comp, prov.all_valid_components)
        print('temp:', temp, 'V/F:', results.vap_frac_mole)

        results = flash_press_prop_2phase(prov, press, 'enthalpy_mole', results.enthalpy_mole, 0, feed_comp, valid=None, previous=None)
        print('temp:', results.temp, 'V/F:', results.vap_frac_mole)
        print('---')

def test_flash_temp_prop():
    components = pure(['ETHANE', 'N-BUTANE', 'N-HEPTANE'])
    prov = PengRobinson(components)
    temp, press = 400, 101325 * 5
    feed_comp = np.array([0.2, 0.3, 0.5])

    results = basic_flash_temp_press_2phase(prov, temp, press, feed_comp, prov.all_valid_components)
    print('press:', press, 'V/F:', results.vap_frac_mole)

    results = flash_temp_prop_2phase(prov, temp, 'enthalpy_mole', results.enthalpy_mole, 0, feed_comp, valid=None, previous=None)
    print('press:', results.press, 'V/F:', results.vap_frac_mole)

def test_flash_prop_vap_frac():
    components = pure(['ETHANE', 'N-BUTANE', 'N-HEPTANE'])
    prov = PengRobinson(components)
    temp, press = 400, 101325 * 10
    feed_comp = np.array([0.2, 0.3, 0.5])

    results = basic_flash_temp_press_2phase(prov, temp, press, feed_comp, prov.all_valid_components)
    print('press:', press, 'temp:', temp, 'V/F:', results.vap_frac_mole)

    results = flash_prop_vap_frac_2phase(prov, 'enthalpy_mole', results.enthalpy_mole, 0, results.vap_frac_mole, feed_comp)
    print('press:', results.press, 'temp:', results.temp, 'V/F:', results.vap_frac_mole)

def test_all():
    press_vs_rho()
    rhos_vs_temp()
    test_flash()
    test_cavett_feed()
    test_flash_press_vap_frac()
    test_flash_temp_vap_frac()
    test_flash_press_prop()
    test_flash_temp_prop()
    test_flash_prop_vap_frac()
    pass

test_all()


text:
#!/bin/env/python

import numpy as np
import matplotlib.pyplot as plt
from sklearn.decomposition import TruncatedSVD

from ..algo.motif import findMotif, findAllMotifInstances
from ..utils.subseq import simMatFromDistTensor

from viz_utils import plotRect, plotRanges

def findAndPlotMotif(seq, lengths, **kwargs):
	motif = findMotif([seq], lengths)
	plotMotif(seq, motif, **kwargs)

def findAndPlotMotifInstances(seq, lengths, truthStartEndPairs=None,
	saveas=None, findMotifKwargs=None, **kwargs):
	# XXX this func will break if seq is a list of seqs, not just one ndarray

	if findMotifKwargs:
		startIdxs, instances, motif = findAllMotifInstances([seq], lengths,
			**findMotifKwargs)
	else:
		startIdxs, instances, motif = findAllMotifInstances([seq], lengths)

	# ------------------------ plot reported motif instances

	endIdxs = startIdxs + motif.length
	startEndPairs = np.c_[startIdxs, endIdxs]
	ax = plotMotifInstances(seq, startEndPairs, **kwargs)

	# ------------------------ plot ground truth
	if truthStartEndPairs is not None and len(truthStartEndPairs):
		try:
			if len(truthStartEndPairs[0]) == 1: # single points, not ranges
				color = 'k'
				truthStartEndPairs = np.asarray(truthStartEndPairs)
				truthStartEndPairs = np.c_[truthStartEndPairs, truthStartEndPairs]
			else:
				color = 'g'
		except: # elements are scalars and so len() throws
			color = 'k'
			truthStartEndPairs = np.asarray(truthStartEndPairs)
			truthStartEndPairs = np.c_[truthStartEndPairs, truthStartEndPairs]

		# make a vert line (spanning less than full graph height)
		# where the labels are
		yMin, yMax = np.min(seq), np.max(seq)
		yRange = yMax - yMin
		lineMin, lineMax = [yMin + frac * yRange for frac in (.4, .6)]
		plotMotifInstances(None, truthStartEndPairs, ax=ax,
			color=color, linestyle='-', lw=2,
			ymin=lineMin, ymax=lineMax)

	if saveas:
		plt.savefig(saveas)
	else:
		plt.show()

def plotMotif(seq, motif, showExtracted=True, color='gray',
	title=None, saveas=None):

	start1 = motif[3]
	start2 = motif[4]
	end1 = start1 + len(motif[0]) - 1
	end2 = start2 + len(motif[1]) - 1

	# just show where the motif is in the original signal
	if not showExtracted:
		_, ax = plt.subplots()
		ax.autoscale(False)
		ax.plot(seq)
		plotRect(ax, start1, end1, color=color)
		plotRect(ax, start2, end2, color=color)
		if saveas:
			plt.savefig(saveas)
		else:
			plt.show()
		return ax

	# set up axes
	ax1 = plt.subplot2grid((2,2), (0,0), colspan=2)
	ax2 = plt.subplot2grid((2,2), (1,0))
	ax3 = plt.subplot2grid((2,2), (1,1))
	ax1.autoscale(tight=True)
	ax2.autoscale(tight=True)
	ax3.autoscale(tight=True)

	# plot raw ts on top and motif instances on the bottom
	ax1.plot(seq, lw=2)
	ax2.plot(motif[0], lw=2)
	ax3.plot(motif[1], lw=2)
	ax1.set_title('Original Signal')
	ax2.set_title('Motif Instance at %d' % start1)
	ax3.set_title('Motif Instance at %d' % start2)

	# draw rects in the ts where the motif is
	plotRect(ax1, start1, end1, color=color)
	plotRect(ax1, start2, end2, color=color)

	plt.tight_layout()
	if saveas:
		plt.savefig(saveas)
	else:
		plt.show()

	return ax1, ax2, ax3

def plotMotifInstances(seq, startEndIdxPairs, title=None, ax=None,
	saveas=None, **kwargs):

	if ax is None:
		_, ax = plt.subplots()
		# ax.autoscale(False) # makes it not actually work...

	if seq is not None and len(seq): # plot original seq if one is provided
		ax.plot(seq, **kwargs)

	plotRanges(ax, startEndIdxPairs, **kwargs)

	if not title:
		title = "Motif Instances in Data"
	ax.set_title(title)
	if seq is not None and len(seq):
		ax.set_ylim([np.min(seq), np.max(seq)])
		ax.set_xlim([0, len(seq)])

	if saveas:
		plt.savefig(saveas)

	return ax

def showPairwiseSims(origSignal, m, simMat, clamp=True, pruneCorrAbove=-1,
	plotMotifs=True, showEigenVect=False, hasPadding=True, saveas=None):

	print "origSignal shape", origSignal.shape
	# padLen = len(origSignal) - simMat.shape[1]
	padLen = m - 1 if hasPadding else 0
	subseqLen = m

	plt.figure(figsize=(8,10))
	if showEigenVect:
		ax1 = plt.subplot2grid((20,10), (0,0), colspan=8, rowspan=5)
		ax2 = plt.subplot2grid((20,10), (5,0), colspan=8, rowspan=15)
		ax3 = plt.subplot2grid((20,10), (5,8), colspan=2, rowspan=15)
		ax3.autoscale(tight=True)
		ax3.set_title('Extracted')
	else:
		ax1 = plt.subplot2grid((4,2), (0,0), colspan=2)
		ax2 = plt.subplot2grid((4,2), (1,0), colspan=2, rowspan=3)
	ax1.autoscale(tight=True)
	ax2.autoscale(tight=True)
	ax1.set_title('Original Signal')
	ax1.set_ylabel('Value')

	if pruneCorrAbove > 0:
		ax2.set_title('Subsequence Cosine Similarities to Dictionary Sequences')
	else:
		ax2.set_title('Subsequence Pairwise Cosine Similarities')
	ax2.set_xlabel('Subsequence Start Index')
	ax2.set_ylabel('"Dictionary" Sequence Number')

	seq = origSignal
	imgMat = simMat

	print "imgMat shape: ", imgMat.shape

	# 	# show magnitude of similarities in each row in descending order; there are
	# 	# only about 60 entries > .01 in *any* row for msrc, and way fewer in most
	# 	# plt.figure()
	# 	# thresh = .5
	# 	# sortedSimsByRow = np.sort(imgMat, axis=1)
	# 	# sortedSimsByRow = sortedSimsByRow[:, ::-1]
	# 	# nonzeroCols = np.sum(sortedSimsByRow, axis=0) > thresh # ignore tiny similarities
	# 	# sortedSimsByRow = sortedSimsByRow[:, nonzeroCols]
	# 	# # plt.imshow(sortedSimsByRow)
	# 	# # plt.plot(np.mean(sortedSimsByRow, axis=1))
	# 	# plt.plot(np.sum(sortedSimsByRow > thresh, axis=1)) # entries > thresh per row

	# if pruneCorrAbove > 0.:
	# 	print "ImgMat Shape:"
	# 	print imgMat.shape
	# 	imgMat = removeCorrelatedRows(imgMat, pruneCorrAbove)
	# 	print imgMat.shape
	# 	print "NaNs at:", np.where(np.isnan(imgMat))[0]
	# 	print "Infs at:", np.where(np.isinf(imgMat))[0]

	# power iteration to see what we get
	if showEigenVect:
		width = int(subseqLen * 1.5)
		nRows, nCols = imgMat.shape
		nPositions = nCols - width + 1
		if nPositions > 1:
			elementsPerPosition = nRows * width # size of 2d slice
			dataMat = np.empty((nPositions, elementsPerPosition))
			# for i in range(nPositions): 		 	# step by 1
			for i in range(0, nPositions, width): 	# step by width, so non-overlapping
				startCol = i
				endCol = startCol + width
				data = imgMat[:, startCol:endCol]
				dataMat[i] = data.flatten()
			# ah; power iteration is for cov matrix, cuz needs a square mat
			# v = np.ones(elementsPerPosition) / elementsPerPosition # uniform start vect
			# for i in range(3):
			# 	v = np.dot(dataMat.T, v)
			svd = TruncatedSVD(n_components=1, random_state=42)
			svd.fit(dataMat)
			v = svd.components_[0]
			learnedFilt = v.reshape((nRows, width))
			ax3.imshow(learnedFilt) # seems to be pretty good
			# plt.show()

	ax1.plot(seq)
	ax2.imshow(imgMat, interpolation='nearest', aspect='auto')
	plt.tight_layout()

	if plotMotifs:
		searchSeq = seq
		print "searchSeq shape:", searchSeq.shape
		motif = findMotif([searchSeq], subseqLen) # motif of min length
		start1 = motif[3]
		start2 = motif[4]
		end1 = start1 + len(motif[0]) - 1
		end2 = start2 + len(motif[1]) - 1

		ax2.autoscale(False)
		color = 'grey'
		plotRect(ax1, start1, end1, color=color)
		plotRect(ax2, start1, end1, color=color)
		plotRect(ax1, start2, end2, color=color)
		plotRect(ax2, start2, end2, color=color)

		print "imgMat shape: ", imgMat.shape
		print "padLen: ", padLen
		if padLen:
			searchSeq = imgMat[:,:-padLen].T
		else:
			searchSeq = imgMat.T
		print "searchSeq shape:", searchSeq.shape
		print "subseqLen:", subseqLen
		motif = findMotif([searchSeq], subseqLen) # motif of min length
		start1 = motif[3]
		start2 = motif[4]
		end1 = start1 + len(motif[0]) - 1
		end2 = start2 + len(motif[1]) - 1
		print [start1, end1, start2, end2]

		color = 'm' # magenta
		plotRect(ax1, start1, end1, color=color)
		plotRect(ax2, start1, end1, color=color)
		plotRect(ax1, start2, end2, color=color)
		plotRect(ax2, start2, end2, color=color)

	if saveas:
		plt.savefig(saveas)
	else:
		plt.show()

	if showEigenVect:
		return ax1, ax2, ax3
	return ax1, ax2

def showPairwiseDists(origSignal, m, Dtensor, **kwargs):
	padLen = len(origSignal) - Dtensor.shape[1]
	simMat = simMatFromDistTensor(Dtensor, m, padLen)
	showPairwiseSims(origSignal, m, simMat, **kwargs)

max_stars_repo_name: icdm-extract/extract

max_issues_repo_name: icdm-extract/extract

max_forks_repo_name: icdm-extract/extract


text:
import requests
import json
import pickle
import numpy as np

class WikipediaRetriever:
    """
    Retreiver using the official Wikipedia API
    """
    def __init__(self, docs, k):
        """
        :param docs: dictionary form of Wikipedia
        :param k:  number of articles to return
        """
        self.S = requests.Session()
        self.URL = "https://ar.wikipedia.org/w/api.php"
        self.k = k
        self.docs = docs

    def get_topk_docs_scores(self, query):
        PARAMS = {
            'action': "query",
            'list': "search",
            'srsearch': query,
            'format': "json"
        }
        R = self.S.get(url=self.URL, params=PARAMS)
        DATA = R.json()
        L = DATA['query']['search']
        article_titles = []
        for i in range(0, len(L)):
            if L[i]['title'] in self.docs:
                article_titles.append(L[i]['title'])
        top_docs = []
        docs_scores = []
        rank = 1
        for title in article_titles[:min(self.k, len(article_titles))]:
            if title in self.docs:
                for par in self.docs[title]:
                    if len(par)>100:
                        top_docs.append(par)
                        docs_scores.append(1/rank)
            rank += 1
        docs_scores = np.asarray(docs_scores)
        docs_scores = docs_scores / docs_scores.sum(axis=0, keepdims=1)
        return top_docs, docs_scores

    def get_topk_docs(self, query):
        PARAMS = {
            'action': "query",
            'list': "search",
            'srsearch': query,
            'format': "json"
        }
        R = self.S.get(url=self.URL, params=PARAMS)
        DATA = R.json()
        L = DATA['query']['search']
        article_titles = [""] * len(L)
        for i in range(0, len(L)):
            article_titles[i] = L[i]['title']
        top_docs = []
        for title in article_titles:
            if title in self.docs:
                top_docs.append(self.docs[title])
        return top_docs[:min(self.k, len(top_docs))]

def test_WikipediaRetriever():
    wiki_data = pickle.load(open("../wiki_extractor/arwiki.p", "rb"))
    r = WikipediaRetriever(wiki_data, 5)
    print(r.get_topk_docs_scores("نادي ليفربول"))

#test_WikipediaRetriever()

max_stars_repo_name: AseelAlshorafa/SOQAL


text:
import os
import pandas as pd
import numpy as np
from datetime import date

class Reimbursements:

    FILE_BASE_NAME = 'reimbursements.xz'

    CSV_PARAMS = {
        'compression': 'xz',
        'encoding': 'utf-8',
        'index': False
    }

    YEARS = [n for n in range(2009, date.today().year+1)]

    def __init__(self, path):
        self.path = path

    def read_csv(self, name):
        filepath = os.path.join(self.path, name)
        print('Loading {}…'.format(name))
        dtype = {
            'applicant_id': np.str,
            'batch_number': np.str,
            'cnpj_cpf': np.str,
            'congressperson_document': np.str,
            'congressperson_id': np.str,
            'document_id': np.str,
            'document_number': np.str,
            'document_type': np.str,
            'leg_of_the_trip': np.str,
            'passenger': np.str,
            'reimbursement_number': np.str,
            'subquota_group_description': np.str,
            'subquota_group_id': np.str,
            'subquota_number': np.str,
            'term_id': np.str,
        }
        return pd.read_csv(filepath, dtype=dtype)

    @property
    def receipts(self):
        print('Merging all datasets…')
        datasets = ["reimbursements-{}.xz".format(n) for n in self.YEARS]
        data = (self.read_csv(name) for name in datasets)
        return pd.concat(data)

    @staticmethod
    def aggregate(grouped, old, new, func):
        """
        Gets a GroupBy object, aggregates it on `old` using `func`, then rename
        the series name from `old` to `new`, returning a DataFrame.
        """
        output = grouped[old].agg(func)
        output = output.rename(index=new, inplace=True)
        return output.reset_index()

    def group(self, receipts):
        print('Dropping rows without document_value or reimbursement_number…')
        subset = ('document_value', 'reimbursement_number')
        receipts = receipts.dropna(subset=subset)

        groupby_keys = ('year', 'applicant_id', 'document_id')
        receipts = receipts.dropna(subset=subset + groupby_keys)

        receipts = receipts[receipts['document_value'] != 0]
        receipts = receipts[receipts['reimbursement_number'] != '0']
        receipts = receipts[receipts['year'] != 0]
        receipts = receipts[receipts['applicant_id'] != '0']
        receipts = receipts[receipts['document_id'] != '0']

        print('Grouping dataset by applicant_id, document_id and year…')
        grouped = receipts.groupby(groupby_keys)

        print('Gathering all reimbursement numbers together…')
        numbers = self.aggregate(
            grouped,
            'reimbursement_number',
            'reimbursement_numbers',
            lambda x: ', '.join(set(x))
        )

        print('Summing all net values together…')
        net_total = self.aggregate(
            grouped,
            'net_value',
            'total_net_value',
            np.sum
        )

        print('Summing all reimbursement values together…')
        total = self.aggregate(
            grouped,
            'reimbursement_value',
            'reimbursement_value_total',
            np.sum
        )

        print('Generating the new dataset…')
        final = pd.merge(
            pd.merge(pd.merge(total, net_total, on=groupby_keys), numbers, on=groupby_keys),
            receipts,
            on=groupby_keys
        )
        final = final.drop_duplicates(subset=groupby_keys)
        final.rename(columns={'net_value': 'net_values',
                              'reimbursement_value': 'reimbursement_values'},
                     inplace=True)
        final = final.drop('reimbursement_number', 1)
        return final

    @staticmethod
    def unique_str(strings):
        return ', '.join(set(strings))

    def write_reimbursement_file(self, receipts):
        print('Casting changes to a new DataFrame…')
        receipts_data_frame = pd.DataFrame(data=receipts)

        print('Writing it to file…')
        filepath = os.path.join(self.path, self.FILE_BASE_NAME)
        receipts_data_frame.to_csv(filepath, **self.CSV_PARAMS)

        print('Done.')

if __name__ == '__main__':
    reimbursements = Reimbursements(os.path.join('tmp', 'serenata-data'))
    reimbursements_group_data_frame = reimbursements.group(reimbursements.receipts)
    reimbursements.write_reimbursement_file(reimbursements_group_data_frame)

max_stars_repo_name: Amandasv/serenata-toolbox

max_issues_repo_name: Amandasv/serenata-toolbox

max_forks_repo_name: Amandasv/serenata-toolbox