text:
from keras.models import Sequential
from keras.layers import Dense, Activation
from keras.optimizers import SGD, Adagrad
from utilityV2.utils import read_data, generate_datasets, evaluate_result
from sklearn.preprocessing import normalize
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np

"""
This version 2 will attempt to use keras to create a multi-layer neural network
to predict the next price of a particular cryptocurrency given the current
price, volume, high price, low price of a given timestep
"""

data = read_data('crypto-test-data-82hrs.csv')
data = data[data['ticker'] == 'ETHBTC']
data['timestamp'] = pd.to_datetime(data.timestamp)
data = data[['high', 'low', 'price', 'volume', 'timestamp']].sort_values(by='timestamp')
data = data.set_index('timestamp')

X_train, X_test, Y_train, Y_test = generate_datasets(data)

print('X_train shape: ' + str(X_train.shape))
print('X_test shape: ' + str(X_test.shape))
print('Y_train shape: ' + str(Y_train.shape))

model = Sequential()
model.add(Dense(16, input_dim=X_train.shape[1], activation='relu'))
model.add(Dense(6, activation='relu'))
model.add(Dense(Y_train.shape[1]))
model.compile(optimizer=SGD(lr=0.03), loss='mse') # Stochastic Gradient Descent
# model.compile(optimizer=Adagrad(lr=0.01), loss='mse') # Adam, to be in sync with v1, use SGD

model.fit(X_train, Y_train, epochs=6000, batch_size=X_train.shape[0]) # To be in sync with v1, no batch training

y_hat_train = model.predict(X_train, batch_size=X_train.shape[0])
evaluate_result(y_hat_train, X_train, Y_train, model, 'train')

y_hat_test = model.predict(X_test, batch_size=X_test.shape[0])
evaluate_result(y_hat_test, X_test, Y_test, model, 'test')

max_stars_repo_name: workofart/ml-trading-playground

max_issues_repo_name: workofart/ml-trading-playground

max_forks_repo_name: workofart/ml-trading-playground


text:
###############################################################################
# WaterTAP Copyright (c) 2021, The Regents of the University of California,
# through Lawrence Berkeley National Laboratory, Oak Ridge National
# Laboratory, National Renewable Energy Laboratory, and National Energy
# Technology Laboratory (subject to receipt of any required approvals from
# the U.S. Dept. of Energy). All rights reserved.
#
# Please see the files COPYRIGHT.md and LICENSE.md for full copyright and license
# information, respectively. These files are also available online at the URL
# "https://github.com/watertap-org/watertap/"
#
###############################################################################

"""Feed blocks for supported property packages"""

from pyomo.environ import Constraint
from idaes.generic_models.unit_models import Feed
from watertap.examples.flowsheets.full_treatment_train.model_components import property_models
from idaes.core.util.scaling import calculate_scaling_factors

def build_feed(m, base='TDS'):
    """
    Build a feed block for a specified property base. The state vars are fixed to the standard condition.
    Property bases include: 'TDS', 'ion', 'salt'
    """

    # feed block
    prop = property_models.get_prop(m, base=base)

    # build
    m.fs.feed = Feed(default={'property_package': prop})

    # specify
    property_models.specify_feed(m.fs.feed.properties[0], base=base)
    m.fs.feed.properties[0].mass_frac_phase_comp  # touch so the block can be initialized

ext: py

max_stars_repo_name: srikanthallu/proteuslib

max_issues_repo_name: srikanthallu/proteuslib

max_forks_repo_name: srikanthallu/proteuslib

reason: from pyomo


text:
import sys
sys.path.insert(0, '..')
import updates_msbm_vi_iter
import updates_msbm_vi
import updates_msbm2_vi_iter
import updates_msbm2_vi

import os
import util
import init_msbm_vi as im
import numpy as np
import numpy.random as npr
import pdb

class TestUpdates:

    def __init__(self, updater_einsum, updater_iter, file_url, remove_self_loops):

        print("\n\n=====================================================================================")
        print("Comparing '{:}' with '{:}' with remove_self_loops={:}"
              .format(updater_einsum.__name__, updater_iter.__name__, remove_self_loops))
        print("=====================================================================================")

        data = util.load_data(file_url)
        self.data = data
        print('')

        K = 8

        data['X'] = data['X'][:K, :, :]
        data['NON_X'] = data['NON_X'][:K, :, :]

        data['K'] = K
        data['Y'] = data['Y'][:K, :]

        prior = dict()
        prior['ALPHA_0'] = 0.5
        prior['BETA_0'] = 0.5
        prior['NU_0'] = 0.5
        prior['ZETA_0'] = 0.5
        self.prior = prior

        # assigning hyper-parameters from ground truth (cheating)
        hyper = dict()
        hyper['M'] = data['M']
        hyper['Q'] = data['Q']
        self.hyper = hyper

        # initialize moments
        mom = dict()
        npr.seed(1)
        mode = 'random'
        mom['ALPHA'] = im.init_ALPHA(data, hyper, mode)
        mom['BETA'] = im.init_BETA(data, hyper, mode)
        mom['NU'] = im.init_NU(data, hyper, mode)
        mom['ZETA'] = im.init_ZETA(data, hyper, mode)
        mom['MU'] = im.init_MU(data, hyper, mode)
        mom['LOG_MU'] = np.log(mom['MU'])
        mom['TAU'] = im.init_TAU(data, hyper, mode)
        mom['LOG_TAU'] = np.log(mom['TAU'])
        self.mom = mom

        par = dict()
        par['MAX_ITER'] = 1000
        par['TOL_ELBO'] = 1.e-16
        par['ALG'] = 'cavi'
        par['kappa'] = 1.0
        self.par = par

        self.remove_self_loops = remove_self_loops
        self.msbm_einsum = updater_einsum
        self.msbm_iter = updater_iter

    def test_update_Pi(self):
        print('--- Pi ---')

        NEW_ALPHA1, NEW_BETA1 = self.msbm_einsum.update_Pi(self.data, self.prior, self.hyper, self.mom, self.par,
                                                         remove_self_loops=self.remove_self_loops)
        NEW_ALPHA2, NEW_BETA2 = self.msbm_iter.update_Pi(self.data, self.prior, self.hyper, self.mom, self.par,
                                                         remove_self_loops=self.remove_self_loops)

        print("--ALPHA:")
        self.eval_diff(NEW_ALPHA1, NEW_ALPHA2)
        print("--BETA:")
        self.eval_diff(NEW_BETA1, NEW_BETA2)

    def test_update_Z(self):
        print('--- Z ---')

        NEW_LOG_TAU1 = self.msbm_einsum.update_Z(self.data, self.prior, self.hyper, self.mom, self.par,
                                                remove_self_loops=self.remove_self_loops)
        NEW_LOG_TAU2 = self.msbm_iter.update_Z(self.data, self.prior, self.hyper, self.mom, self.par,
                                               remove_self_loops=self.remove_self_loops)

        self.eval_diff(NEW_LOG_TAU1, NEW_LOG_TAU2)

    def test_update_Y(self):
        print('--- Y ---')

        NEW_LOG_MU1 = self.msbm_einsum.update_Y(self.data, self.prior, self.hyper, self.mom, self.par,
                                                remove_self_loops=self.remove_self_loops)
        NEW_LOG_MU2 = self.msbm_iter.update_Y(self.data, self.prior, self.hyper, self.mom, self.par,
                                              remove_self_loops=self.remove_self_loops)

        self.eval_diff(NEW_LOG_MU1, NEW_LOG_MU2)

    def test_update_gamma(self):
        print('---  Gamma ---')

        NEW_NU1 = self.msbm_einsum.update_gamma(self.data, self.prior, self.hyper, self.mom, self.par,
                                                remove_self_loops=self.remove_self_loops)
        NEW_NU2 = self.msbm_iter.update_gamma(self.data, self.prior, self.hyper, self.mom, self.par)

        self.eval_diff(NEW_NU1, NEW_NU2)

    def test_update_rho(self):
        print('---  Rho ---')

        NEW_ZETA1 = self.msbm_einsum.update_rho(self.data, self.prior, self.hyper, self.mom, self.par,
                                                remove_self_loops=self.remove_self_loops)
        NEW_ZETA2 = self.msbm_iter.update_rho(self.data, self.prior, self.hyper, self.mom, self.par)

        self.eval_diff(NEW_ZETA1, NEW_ZETA2)

    def test_all(self):
        self.test_update_Pi()
        self.test_update_Z()
        self.test_update_Y()
        self.test_update_gamma()
        self.test_update_rho()

    def eval_diff(self, X1, X2):

        ind = np.unravel_index(np.argmax(X1 - X2), X1.shape)
        diff = (X1 - X2).ravel()
        print('Mean abs entry-wise error: {:1.3e}'.format(np.mean(np.abs(diff))))
        print('Max abs entry-wise error:  {:1.3e}'.format(np.max(np.abs(diff))))
        print('On entry with index:')
        print(ind)

# ###########################################################
# ###########################################################
# ###########################################################

if __name__ == '__main__':

    file_url = os.path.join('..', 'experiments', 'two_prototype', 'data', 'twoprototype_105_250.pickle')

    remove_self_loops = False

    updater_einsum = updates_msbm_vi
    updater_iter = updates_msbm_vi_iter

    runner = TestUpdates(updater_einsum, updater_iter, file_url, remove_self_loops)
    runner.test_all()

    updater_einsum = updates_msbm2_vi
    updater_iter = updates_msbm2_vi_iter

    runner = TestUpdates(updater_einsum, updater_iter, file_url, remove_self_loops)
    runner.test_all()

    remove_self_loops = True

    updater_einsum = updates_msbm_vi
    updater_iter = updates_msbm_vi_iter

    runner = TestUpdates(updater_einsum, updater_iter, file_url, remove_self_loops)
    runner.test_all()

    updater_einsum = updates_msbm2_vi
    updater_iter = updates_msbm2_vi_iter

    runner = TestUpdates(updater_einsum, updater_iter, file_url, remove_self_loops)
    runner.test_all()

max_stars_repo_name: cassianobecker/msbm

max_issues_repo_name: cassianobecker/msbm

max_forks_repo_name: cassianobecker/msbm

reason: import numpy


text:
import dash
import dash_html_components as html
import dash_core_components as dcc
from dash.dependencies import Input, Output

import pandas as pd
import numpy as np

from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
import sklearn.datasets

import plotly.graph_objs as go

wine_data = sklearn.datasets.load_wine()
Y = wine_data["target"]
wine_data = pd.DataFrame(wine_data["data"], columns=wine_data["feature_names"])
wine_data["target"] = Y
scaler = StandardScaler()
wine_scaled = scaler.fit_transform(wine_data.drop(["target"], axis = 1))
wine_scaled = pd.DataFrame(wine_scaled, columns=wine_data.drop(["target"], axis = 1).columns)
pca = PCA(random_state=42)
principalComponents = pca.fit_transform(wine_scaled)
comp_df = pd.DataFrame(principalComponents)
comp_df["target"] = wine_data["target"]

def get_plot(num_value):
	if num_value == 3:
		data = []
		for i in [0, 1, 2]:
		    trace = go.Scatter3d(x = comp_df[comp_df["target"] == i][0],
		                         y = comp_df[comp_df["target"] == i][1],
		                         z = comp_df[comp_df["target"] == i][2],
		                         mode = 'markers',
		                         name = "class_{}".format(i),
		                         marker = dict(size = 3))
		    data.append(trace)

		figure = dict(data=data)
	else:
		data = []
		for i in [0, 1, 2]:
		    trace = go.Scatter(x = comp_df[comp_df["target"] == i][0],
		                         y = comp_df[comp_df["target"] == i][1],
		                         mode = 'markers',
		                         name = "class_{}".format(i),
		                         marker = dict(size = 3))
		    data.append(trace)

		figure = dict(data=data)
	return figure

figure_ = get_plot(num_value = 3)

exp_var = pca.explained_variance_ratio_
exp_var_trace = go.Line(y = exp_var.cumsum(),
               x = np.arange(1, len(exp_var) + 1))
exp_var_data = [exp_var_trace]
exp_var_figure = dict(data = exp_var_data)

app = dash.Dash()
app.css.append_css({'external_url': 'https://codepen.io/chriddyp/pen/bWLwgP.css'})

app.layout = html.Div([
	html.Div([dcc.Dropdown(
        id = 'graph_option',
        options=[
            {'label': '3D', 'value': 3},
            {'label': '2D', 'value': 2},
        ],
        value = 3)], className = "row"),

			html.Div([
					html.Div([dcc.Graph(id='exp_var', figure = exp_var_figure)], className = "six columns"),
					html.Div([dcc.Graph(id='scatter', figure = figure_)], className = "six columns")
				], className = "row")],

	className = "container")

@app.callback(
    Output(component_id='scatter', component_property='figure'),
    [Input(component_id='graph_option', component_property='value')]
)
def update_text_output(value_):
    return get_plot(num_value = value_)

if __name__ == '__main__':
    app.run_server(debug=True)


text:
from __future__ import print_function

import numpy as np
from openmdao.core.component import Component
from openmdao.api import IndepVarComp, Component, Problem, Group
class PodGeometry(Component):
    """
    Notes
    ------
    Computes to corss sectional area, length, and planform area of the pod based on the sizes of internal components
    and the necessary duct area within the pod based on compressor peformance.  Assumes isentropic compression and a
    compressor exit mach number of .3. Also calculate blockage factors based on pressurized cylinder equations.

    Params
    ------
    A_payload : float
        Cross sectional area of passenger compartment. Default value is 2.72
    gam : float
        Ratio of specific heats. Default value is 1.4
    R : float
        Ideal gas constant. Default valut is 287 J/(m*K).
    p_tunnel : float
        Pressure of air in tube.  Default value is 850 Pa.  Value will come from vacuum component
    M_pod : float
        pod Mach number. Default value is .8. value will be set by user
    A_p : float
        Cross sectional area of passenger compartment. Default value is 2.72 m**2
    L_comp : float
        length of the compressor. Default value is 1.0 m.
    L_bat : float
        length of battery. Default value is 1.0 m.
    L_inverter : float
        length of inverter. Default value is 1.0 m.
    L_trans : float
        length of transformer. Default value is 1.0 m
    L_p : float
        length of passenger compartment. Default value is 11.2 m
    L_conv : float
        length of the converging section of the nozzle. Default value is .3 m
    L_div : float
        length of the diverging section of the nozzle. Default value is 1.5 m
    L_motor : float
        length of the motor. Default value is 1.0 m
    p_duct : float
        Static pressure in the duct. Default value is 6800.0 Pa
    p_passenger : float
        Static pressure in passenger section. Default value is 101.3e3 Pa
    rho_pod : float
        Density of pod material. Default value is 2700.0 kg/m**3
    n_passengers : float
        Number of passengers per pod. Default value is 28
    SF : float
        Structural safety factor. Default value is 1.5
    Su : float
        Ultimate strength of pod material. Default value is 50.0e6 Pa
    A_duct : float
        Area of flow duct within pod. Default value is .3 m**2
    dl_passenger : float
        Length of passenger compartment per passenger row. Default value is .8 m.
    g : float
        gravitational acceleration. Default value is 9.81 m/s**2

    Returns
    -------
    A_pod : float
        Cross sectional area of pod.
    D_pod : float
        Diameter of pod
    S : float
        Platform area of the pod
    L_pod : float
        Length of Pod
    t_passenger : float
        thickness of structure in passenger section.
    t_pod : float
        thickness of outer pod section.
    BF : float
        Pod blockage factor.
    beta : float
        Duct blockage factor.
    """

    def __init__(self):
        super(PodGeometry, self).__init__()

        self.add_param('L_comp', val=1.0, desc='Length of Compressor', units='m')
        self.add_param('L_bat', val=1.0, desc='Length of Battery', units='m')
        self.add_param('L_motor', val=1.0, desc='Length of Motor', units='m')
        self.add_param('L_inverter', val=0.0, desc='Length of Inverter', units='m')
        self.add_param('L_trans', val=1.0, desc='Length of Transformer', units='m')
        self.add_param('L_p', val=11.2, desc='Payload Length', units='m')
        self.add_param('L_conv', val=.3, desc='Converging Lenth', units='m')
        self.add_param('L_div', val=1.5, desc='Diverging Length', units='m')
        self.add_param('L_inlet', 2.5, desc = 'Inlet length', units = 'm')
        self.add_param('p_tunnel', val=850.0, desc='tunnel pressure', units='Pa')
        self.add_param('A_payload', val=2.72, desc='Cross sectional area of passenger compartment')
        self.add_param('p_duct', 6800.0, desc = 'duct pressure', units = 'Pa')
        self.add_param('p_passenger', 101.0e3, desc = 'Passenger compartment pressure', units = 'Pa')
        self.add_param('rho_pod', 2700.0, desc = 'material density', units = 'kg/m**3')
        self.add_param('n_passengers', 28., desc = 'number of passengers', units = 'unitless')
        self.add_param('dm_passenger', 166.0, desc = 'mass per passenger', units = 'kg')
        self.add_param('SF', 1.5, desc = 'safety factor for pressure cylinder', units = 'unitless')
        self.add_param('Su', 50.0e6, desc = 'ultimate strength', units = 'Pa')
        self.add_param('A_duct', .3, desc = 'tunnel pressure', units = 'm**2')
        self.add_param('dl_passenger', .8, desc = 'passenger compartment length per person', units = 'm')
        self.add_param('g', 9.81, desc = 'gravity', units = 'm/s**2')

        self.add_output('A_pod', val = 0.0, desc = 'cross sectional area of pod', units = 'm**2')
        self.add_output('D_pod', val=0.0, desc='pod diameter', units='m')
        self.add_output('S', val=0.0, desc='platform area of pod', units='m**2')
        self.add_output('L_pod', val=0.0, desc='length of pod', units='m')
        self.add_output('t_passenger', 0.0, desc = 'passenger compartment thickness', units = 'm')
        self.add_output('t_pod', 0.0, desc = 'outer pod thickness', units = 'm')
        self.add_output('BF', 0.0, desc = 'Tunnel blockage Factor', units = 'unitless')
        self.add_output('beta', 0.0, desc = 'duct blockage factor', units = 'unitless')

    def solve_nonlinear(self, p, u, r):

        dp_passenger = p['p_passenger'] - p['p_duct']                   #Calculate pressure differential across passenger section wall
        r_passenger = np.sqrt(p['A_payload']/np.pi)                     #Calculate radius of the passenger section
        t_passenger = (dp_passenger*r_passenger)/(p['Su']/p['SF'])      #Calculate pod thickness based on pressurized cylinder equations

        L_passenger = (p['n_passengers']*p['dl_passenger'])/2.0         #Calculate length of passenger section
        r_pod = np.sqrt((p['A_duct']+p['A_payload'])/np.pi)             #Calculate pod radius initial value for sizing

        #Calculate mass of passenger section. Factor of 1.5 is to account for structure, luggage, etc.
        m_passenger = (p['n_passengers']*p['dm_passenger']*1.5) + p['rho_pod']*np.pi*(((r_passenger+t_passenger)**2)-(r_passenger**2))*L_passenger

        dx = (m_passenger*p['g'])/((p['Su']/5.0)*L_passenger)           #Calculate thickness of structure supporting passenger section according to yield condition with safety factor = 5
        A_cross = dx*(r_pod-r_passenger)                                #Calculate cross sectional area of passenger compartment support structure

        beta = (A_cross + np.pi*(((r_passenger+t_passenger)**2) - (r_passenger**2)))/p['A_payload']     #Calculate duct blockage factor

        dp_pod = p['p_passenger'] - p['p_tunnel']                       #Calculate pressure differntial assuming atomospheric internal pressure
        t_pod = (dp_pod*r_pod)/(p['Su']/p['SF'])                        #Calculate pod thickness based on pressurized cylinder equationss
        r_pod = np.sqrt((p['A_duct']+(1.0+beta)*p['A_payload'])/np.pi)  #Calculate corrected value of pod radius to account for beta

        BF = (p['A_payload']*(1.0+beta) + p['A_duct'])/(np.pi*((r_pod+t_pod)**2.0))                     #Calculate pod blockage factor
        A_pod = (p['A_duct'] + (1+beta)*p['A_payload'])/BF              #Calculate cross sectional area of the pod
        D_pod = np.sqrt((4*A_pod)/np.pi)                                #Calculate pod diameter

        L_pod = p['L_inlet'] + p['L_comp'] + p['L_bat'] + p['L_motor'] + p['L_inverter'] + p['L_trans'] + p['L_p'] + p['L_conv'] + p['L_div'] + L_passenger  #Calculate pod length
        S = D_pod*L_pod                                                 #Calculate pod planform area

        u['A_pod'] = A_pod
        u['D_pod'] = D_pod
        u['L_pod'] = L_pod
        u['S'] = S
        u['t_passenger'] = t_passenger
        u['t_pod'] = t_pod
        u['BF'] = BF
        u['beta'] = beta

if __name__ == '__main__':
    top = Problem()
    root = top.root = Group()

    root.add('p', PodGeometry())

    top.setup()
    top.run()

    print('\n')
    print('Pod cross section = %f m^2' % top['p.A_pod'])
    print('Pod Diameter = %f m' % top['p.D_pod'])
    print('Pod Length = %f m' % top['p.L_pod'])
    print('Pod planform area = %f m^2' % top['p.S'])
    print('Passenger compartment thickness = %f m' % top['p.t_passenger'])
    print('Pod thickness is = %f m' % top['p.t_pod'])
    print('Pod blockage factor = %f' % top['p.BF'])
    print('Duct blockage factor is = %f' % top['p.beta'])


text:
from __future__ import division, print_function
import numpy as np
from scipy.stats import gamma, exponnorm

class synthetic_light_curve_generator:
    def __init__(self, T, N, rseed=None):
        """
        Class for simple synthetic light curve generation

        Light curves are time series of stellar magnitude or flux as a
        function of time. In the following we consider Earth-based surveys.
        Light curves are irregularly sampled because the measurements are
        not taken at the same time every night due to observation
        constraints. Light curves may have data gaps, i.e. periods of time
        where no observations were registered. The noise present in the
        light curves comes from the sky background, the Earth atmosphere,
        telescope systematics, etc.

        Parameters
        ----------
        T: float
            Time span of the vector, i.e. how long it is in time
        N: positive integer
            Number of samples of the resulting time vector
        rseed:
            Random seed to feed the random number generator

        """
        self.rseed = rseed
        self.t = irregular_sampling(T, N, rseed)
        self.A = 1.0

    def set_model(self, f0, A, model='trigonometric'):
        """
        See also
        --------
        trigonometric_model

        """
        self.f0 = f0
        self.y_clean, vary = trigonometric_model(self.t, f0, A)
        # y_clean is normalized to have unit standard deviation
        self.y_clean = self.y_clean/np.sqrt(vary)

    def get_fundamental_frequency(self):
        return self.f0

    def get_clean_signal(self):
        return self.t, self.A*self.y_clean

    def draw_noisy_time_series(self, SNR=1.0, red_noise_ratio=0.25, outlier_ratio=0.0):
        """
        A function to draw a noisy time series based on the clean model
        such that
            y_noisy = y + yw + yr,
        where yw is white noise, yr is red noise and y will be rescaled
        so that y_noisy complies with the specified signal-to-noise
        ratio (SNR).

        Parameters
        ---------
        SNR: float
            Signal-to-noise ratio of the resulting contaminated signal in
            decibels [dB]. SNR is defined as
                SNR = 10*log(var_signal/var_noise), hence
            NR var_signal/var_noise
                 10  10
                  7   5
                  3   2
                  0   1
                 -3   0.5
                 -7   0.2
                -10   0.1
        red_noise_variance: float in [0, 1]
            The variance of the red noise component is set according to
            Var(yw)*red_noise_ratio. Set this to zero to obtain uncertainties
            that explain the noise perfectly
        outlier_ratio: float in [0, 1]
            Percentage of outlier data points

        Returns
        -------
        t: ndarray
            Vector containing the time instants
        y_noisy: ndarray
            Vector containing the contaminated signal
        s: ndarray
            Vector containing the uncertainties associated to the white
            noise component

        """
        if outlier_ratio < 0.0 or outlier_ratio > 1.0:
            raise ValueError("Outlier ratio must be in [0, 1]")
        if red_noise_ratio < 0.0:
            raise ValueError("Red noise ratio must be positive")
        np.random.seed(self.rseed)
        t = self.t
        y_clean = self.y_clean
        N = len(t)
        # First we generate s
        s, mean_s_squared = generate_uncertainties(N, rseed=self.rseed)
        #print(mean_s_squared)
        #print(np.mean(s**2))
        # Draw a heteroscedastic white noise vector
        white_noise = np.random.multivariate_normal(np.zeros(N,), np.diag(s**2))
        # Now we generate a colored noise vector which is unaccounted by s
        red_noise_variance = mean_s_squared*red_noise_ratio
        # First order markovian process to generate
        red_noise = first_order_markov_process(t, red_noise_variance, 1.0, rseed=self.rseed)

        # The following is not ok for irregularly sampled time series because
        # it assumes constant dt=1
        #phi=0.5
        #red_noise = np.random.randn(N)*np.sqrt(red_noise_variance)
        #for i in range(1, N):
        #    red_noise[i] = phi*red_noise[i-1] + np.sqrt(1 - phi**2)*red_noise[i]

        # The final noise vector
        #print("%f %f" % (np.var(white_noise)*red_noise_ratio, np.var(red_noise)))
        noise = white_noise + red_noise
        var_noise = mean_s_squared + red_noise_variance
        SNR_unitless = 10.0**(SNR/10.0)
        self.A = np.sqrt(SNR_unitless*var_noise)
        y = self.A*y_clean
        y_noisy = y + noise
        # Add outliers with a certain percentage
        rperm = np.where(np.random.uniform(size=N) < outlier_ratio)[0]
        outlier = np.random.uniform(5.0*np.std(y), 10.0*np.std(y), size=len(rperm))
        y_noisy[rperm] += outlier
        return t, y_noisy, s

def irregular_sampling(T, N, rseed=None):
    """
    Generates an irregularly sampled time vector by perturbating a
    linearly spaced vector and latter deleting a certain number of
    points

    Parameters
    ----------
    T: float
        Time span of the vector, i.e. how long it is in time
    N: positive integer
        Number of samples of the resulting time vector
    rseed:
        Random seed to feed the random number generator

    Returns
    -------
    t_irr: ndarray
        An irregulary sampled time vector

    """
    sampling_period = (T/float(N))
    N = int(N)
    np.random.seed(rseed)
    t = np.linspace(0, T, num=5*N)
    # First we add jitter
    t[1:-1] += sampling_period*0.5*np.random.randn(5*N-2)
    # Then we do a random permutation and keep only N points
    P = np.random.permutation(5*N)
    t_irr = np.sort(t[P[:N]])
    return t_irr

def trigonometric_model(t, f0, A):
    """
    Generates a simple trigonometric model based on a sum of sine waves

    Parameters
    ----------
    t: ndarray
        Vector containing the time instants where the model will be sampled
    f0: float
        Fundamental frequency of the model
    A: ndarray
        Vector containing the amplitudes of the harmonics in the model.
        The lenght of A defines the number of harmonics

    Returns
    -------
    y: ndarray
        Trigonometric model based
    var_y: float
        Analytical variance of the trigonometric model

    Example
    -------
    Assumming a time vector t has already been specified

    >>> import P4J
    >>> y_clean = P4J.trigonometric_model(t, f0=2.0, A=[1.0, 0.5, .25])

    """
    y = 0.0
    M = len(A)
    var_y = 0.0
    for k in range(0, M):
        y += A[k]*np.sin(2.0*np.pi*t*f0*(k+1))
        var_y += 0.5*A[k]**2
    return y, var_y

def first_order_markov_process(t, variance, time_scale, rseed=None):
    """
    Generates a correlated noise vector using a multivariate normal
    random number generator with zero mean and covariance

        Sigma_ij = s^2 exp(-|t_i - t_j|/l),

    where s is the variance and l is the time scale.
    The Power spectral density associated to  this covariance is

        S(f) = 2*l*s^2/(4*pi^2*f^2*l^2 +1),

    red noise spectrum is defined as proportional to 1/f^2. This
    covariance is the one expected from a first order markov process
    (Reference?)

    Parameters
    ---------
    t: ndarray
        A time vector for which the red noise vector will be sampled
    variance: positive float
        variance of the resulting red noise vector
    time_scale: positive float
        Parameter of the covariance matrix

    Returns
    -------
    red_noise: ndarray
        Vector containing the red noise realizations

    See also
    --------
    power_law_noise

    """
    if variance < 0.0:
        raise ValueError("Variance must be positive")
    if time_scale < 0.0:
        raise ValueError("Time scale must be positive")
    np.random.seed(rseed)
    N = len(t)
    mu = np.zeros(shape=(N,))
    if variance == 0.0:
        return mu
    dt = np.repeat(np.reshape(t, (1, -1)), N, axis=0)
    dt = np.absolute(dt - dt.T)  # This is NxN
    S = variance*np.exp(-np.absolute(dt)/time_scale)
    red_noise = np.random.multivariate_normal(mu, S)
    return red_noise

def power_law_noise(t, variance):
    """
    CAUTION: NOT THOROUGHLY TESTED
    Generates a red noise vector by first generating a covariance
    function based on the expected red noise spectra. Red noise spectrum
    follows a power law

        S(f) = A*f^(-2),

    where A is a constant. Pink noise can be obtained if the exponent is
    changed to -1, and white noise if it is changed to 0.

    This Power Spectral density (PSD) is real and even then

        r(tau) = AT sum cos(2pi k tau/T)/k**2,

    where we assume f[k] = k*df, df = 1/T, Fs = N/T.

    The basel problem gives us: sum (1/k**2) approx pi**2/6, hence we can
    set c according to the desired variance at lag=0

    After the covariance is computed we can draw from a multivariate
    normal distribution to obtain a noise vector

    Parameters
    ---------
    t: ndarray
        A time vector for which the red noise vector will be sampled
    variance: positive float
        variance of the resulting red noise vector

    Returns
    -------
    red_noise: ndarray
        Vector containing the red noise realizations

    See also
    --------
    first_order_markov_process

    """
    if variance < 0.0:
        raise ValueError("Variance must be positive")
    mu = np.zeros(shape=(N,))
    if variance == 0.0:
        return mu
    N = len(t)
    T = (t[-1] - t[0])
    c = (6.0*variance)/(T*np.pi**2)
    f = np.arange(1.0/T, 0.5*N/T, step=1.0/T)  # We ommit f=0.0
    k = f*T
    dt = np.repeat(np.reshape(t, (1, -1)), N, axis=0)
    dt = np.absolute(dt - dt.T)  # This is NxN
    S = np.zeros(shape=(N, N))
    for i in range(0, N):
        for j in range(i, N):
            S[i, j] = np.sum(np.cos(2.0*np.pi*k*dt[i, j]/T)/k**2)*T*c
            S[j, i] = S[i, j]
    red_noise = np.random.multivariate_normal(mu, S)
    return red_noise

def generate_uncertainties(N, dist='Gamma', rseed=None):
    """
    This function generates a uncertainties for the white noise component
    in the synthetic light curve.

    Parameters
    ---------
    N: positive integer
        Lenght of the returned uncertainty vector
    dist: {'EMG', 'Gamma'}
        Probability density function (PDF) used to generate the
        uncertainties
    rseed:
        Seed for the random number generator

    Returns
    -------
    s: ndarray
        Vector containing the uncertainties
    expected_s_2: float
        Expectation of the square of s computed analytically

    """
    np.random.seed(rseed)
    #print(dist)
    if dist == 'EMG':  # Exponential modified Gaussian
        # the mean of a EMG rv is mu + 1/(K*sigma)
        # the variance of a EMG rv is sigma**2 + 1/(K*sigma)**2
        K = 1.824328605481941
        sigma = 0.05*0.068768312946785953
        mu = 0.05*0.87452567616276777
        # IMPORTANT NOTE
        # These parameters were obtained after fitting uncertainties
        # coming from 10,000 light curves of the VVV survey
        expected_s_2 = sigma**2 + mu**2 + 2*K*mu*sigma + 2*K**2*sigma**2
        s = exponnorm.rvs(K, loc=mu, scale=sigma, size=N)
    elif dist == 'Gamma':
        # The mean of a gamma rv is k*sigma
        # The variance of a gamma rv is k*sigma**2
        k = 3.0
        sigma = 0.05/k  #  mean=0.05, var=0.05**2/k
        s = gamma.rvs(k, loc=0.0, scale=sigma, size=N)
        expected_s_2 = k*(1+k)*sigma**2
    return s, expected_s_2

reason: import numpy,from scipy


text:
import numpy as np
from kaldo.observables.observable import Observable

class PhysicalMode(Observable):

    def __init__(self, frequency, min_frequency=0, max_frequency=None, is_nw=False):
        self.frequency = frequency
        self.n_phonons = self.frequency.size
        self.min_frequency = min_frequency
        self.max_frequency = max_frequency
        self.is_nw = is_nw

    def calculate(self, is_loading=False):
        physical_mode = np.ones_like(self.frequency.reshape(self.n_phonons), dtype=bool)
        if self.min_frequency != 0:
            physical_mode = physical_mode & (self.frequency.reshape(self.n_phonons) > self.min_frequency)
        if self.max_frequency is not None:
            physical_mode = physical_mode & (self.frequency.reshape(self.n_phonons) < self.max_frequency)
        if self.is_nw:
            physical_mode[:4] = False
        else:
            physical_mode[:3] = False
        return physical_mode

max_stars_repo_name: kcbhamu/kaldo

max_issues_repo_name: kaituohuo/kaldo

max_forks_repo_name: kaituohuo/kaldo


text:
#!/usr/bin/env python
"""Quantum field theory example

* http://en.wikipedia.org/wiki/Quantum_field_theory

This particular example is a work in progress. Currently it calculates the
scattering amplitude of the process:

    electron + positron -> photon -> electron + positron

in QED (http://en.wikipedia.org/wiki/Quantum_electrodynamics). The aim is to be
able to do any kind of calculations in QED or standard model in SymPy, but
that's a long journey.
"""

from sympy import Basic,exp,Symbol,sin,Rational,I,Mul, Matrix, \
    ones, sqrt, pprint, simplify, Eq, sympify

from sympy.physics import msigma, mgamma

#gamma^mu
gamma0=mgamma(0)
gamma1=mgamma(1)
gamma2=mgamma(2)
gamma3=mgamma(3)
gamma5=mgamma(5)

#sigma_i
sigma1=msigma(1)
sigma2=msigma(2)
sigma3=msigma(3)

E = Symbol("E", real=True)
m = Symbol("m", real=True)

def u(p,r):
    """ p = (p1, p2, p3); r = 0,1 """
    assert r in [1,2]
    p1,p2,p3 = p
    if r == 1:
        ksi = Matrix([ [1],[0] ])
    else:
        ksi = Matrix([ [0],[1] ])
    a = (sigma1*p1 + sigma2*p2 + sigma3*p3) / (E+m) * ksi
    if a ==0:
        a = zeros((2, 1))
    return sqrt(E+m) * Matrix([ [ksi[0,0]], [ksi[1,0]], [a[0,0]], [a[1,0]] ])

def v(p,r):
    """ p = (p1, p2, p3); r = 0,1 """
    assert r in [1,2]
    p1,p2,p3 = p
    if r == 1:
        ksi = Matrix([ [1],[0] ])
    else:
        ksi = -Matrix([ [0],[1] ])
    a = (sigma1*p1 + sigma2*p2 + sigma3*p3) / (E+m) * ksi
    if a ==0:
        a = zeros((2,1))
    return sqrt(E+m) * Matrix([ [a[0,0]], [a[1,0]], [ksi[0,0]], [ksi[1,0]] ])

def pslash(p):
    p1,p2,p3 = p
    p0 = sqrt(m**2+p1**2+p2**2+p3**2)
    return gamma0*p0-gamma1*p1-gamma2*p2-gamma3*p3

def Tr(M):
    return M.trace()

def xprint(lhs, rhs):
    pprint( Eq(sympify(lhs), rhs ) )

def main():
    a=Symbol("a", real=True)
    b=Symbol("b", real=True)
    c=Symbol("c", real=True)

    p = (a,b,c)

    assert u(p, 1).D * u(p, 2) == Matrix(1, 1, [0])
    assert u(p, 2).D * u(p, 1) == Matrix(1, 1, [0])

    p1,p2,p3 =[Symbol(x, real=True) for x in ["p1","p2","p3"]]
    pp1,pp2,pp3 =[Symbol(x, real=True) for x in ["pp1","pp2","pp3"]]
    k1,k2,k3 =[Symbol(x, real=True) for x in ["k1","k2","k3"]]
    kp1,kp2,kp3 =[Symbol(x, real=True) for x in ["kp1","kp2","kp3"]]

    p = (p1,p2,p3)
    pp = (pp1,pp2,pp3)

    k = (k1,k2,k3)
    kp = (kp1,kp2,kp3)

    mu = Symbol("mu")

    e = (pslash(p)+m*ones(4))*(pslash(k)-m*ones(4))
    f = pslash(p)+m*ones(4)
    g = pslash(p)-m*ones(4)

    #pprint(e)
    xprint( 'Tr(f*g)', Tr(f*g) )
    #print Tr(pslash(p) * pslash(k)).expand()

    M0 = [ ( v(pp, 1).D * mgamma(mu) * u(p, 1) ) * ( u(k, 1).D * mgamma(mu,True) * \
             v(kp, 1) ) for mu in range(4)]
    M = M0[0]+M0[1]+M0[2]+M0[3]
    M = M[0]
    assert isinstance(M, Basic)
    #print M
    #print simplify(M)

    d=Symbol("d", real=True) #d=E+m

    xprint('M', M)
    print "-"*40
    M = ((M.subs(E,d-m)).expand() * d**2 ).expand()
    xprint('M2', 1/(E+m)**2 * M)
    print "-"*40
    x,y= M.as_real_imag()
    xprint('Re(M)', x)
    xprint('Im(M)', y)
    e = x**2+y**2
    xprint('abs(M)**2', e)
    print "-"*40
    xprint('Expand(abs(M)**2)', e.expand())

    #print Pauli(1)*Pauli(1)
    #print Pauli(1)**2
    #print Pauli(1)*2*Pauli(1)

if __name__ == "__main__":
    main()

max_stars_repo_name: pernici/sympy

max_issues_repo_name: jegerjensen/sympy

max_forks_repo_name: jegerjensen/sympy

reason: from sympy


text:
from typing import Type, Dict

import numpy as np
from numpy import ndarray

from ..mesh import Mesh, MeshType

class MeshLine(Mesh):
    """One-dimensional mesh."""

    refdom: str = "line"
    brefdom: str = "point"
    meshio_type: str = "line"
    name: str = "One-dimensional"

    p = np.zeros((1, 0), dtype=np.float64)
    t = np.zeros((2, 0), dtype=np.int64)
    t2f = np.zeros((1, 0), dtype=np.int64)
    f2t = np.zeros((2, 0), dtype=np.int64)

    def __init__(self,
                 p: ndarray = None,
                 t: ndarray = None,
                 boundaries: Dict[str, ndarray] = None,
                 subdomains: Dict[str, ndarray] = None,
                 validate: bool = True):
        """Initialize one-dimensional mesh.

        If no parameters are given, initializes a mesh with the nodes 0 and 1.

        Parameters
        ----------
        p
            The nodes of the mesh.
        t
            Optional element connectivity. By default, this is constructed
            using the assumption that `p` is sorted.
        subdomains
            Named subsets of elements.
        boundaries
            Named subsets of boundary facets.
        validate
            If `True`, perform mesh validity checks.

        """
        if p is None and t is None:
            p = np.array([[0., 1.]], dtype=np.float_)
        if len(p.shape) == 1:
            p = np.array([p])
        self.p = p
        self.boundaries = boundaries
        self.subdomains = subdomains

        self.facets = np.arange(self.p.shape[1])[None, :]
        self.t = np.vstack([self.facets[0, :-1],
                            self.facets[0, 1:]]) if t is None else t
        super(MeshLine, self).__init__()
        self._build_mappings()

        if validate:
            self._validate()

    @classmethod
    def init_refdom(cls: Type[MeshType]) -> MeshType:
        """Initialise a mesh consisting of the reference interval [0,1]."""
        return cls()

    def _build_mappings(self):
        """Build t2f and f2t"""

        self.t2f = self.t
        # build f2t
        e_tmp = np.hstack((self.t2f[0, :], self.t2f[1, :]))
        t_tmp = np.tile(np.arange(self.t.shape[1]), (1, 2))[0]

        e_first, ix_first = np.unique(e_tmp, return_index=True)
        e_last, ix_last = np.unique(e_tmp[::-1], return_index=True)
        ix_last = e_tmp.shape[0] - ix_last - 1

        self.f2t = np.zeros((2, self.facets.shape[1]), dtype=np.int64)
        self.f2t[0, e_first] = t_tmp[ix_first]
        self.f2t[1, e_last] = t_tmp[ix_last]

        # second row to zero if repeated (i.e., on boundary)
        self.f2t[1, np.nonzero(self.f2t[0, :] == self.f2t[1, :])[0]] = -1

    def _adaptive_refine(self, marked):
        """Perform an adaptive refine which splits each marked element into
        two."""

        t = self.t
        p = self.p

        mid = range(len(marked)) + np.max(t) + 1

        nonmarked = np.setdiff1d(np.arange(t.shape[1]), marked)

        newp = np.hstack((p, 0.5*(p[:, self.t[0, marked]] +
                                  p[:, self.t[1, marked]])))
        newt = np.vstack((t[0, marked], mid))
        newt = np.hstack((t[:, nonmarked],
                          newt,
                          np.vstack((mid, t[1, marked]))))

        # update fields
        self.p = newp
        self.t = newt

        self._build_mappings()

    def nodes_satisfying(self, test):
        """Return nodes that satisfy some condition.

        Parameters
        ----------
        test : lambda function (1 param)
            Evaluates to 1 or True for nodes belonging
            to the output set.

        """
        return np.nonzero(test(self.p[0, :]))[0]

    def _uniform_refine(self):
        """Perform a single mesh refine that halves 'h'."""
        # rename variables
        p = self.p

        # new vertices and elements
        newp = np.hstack((p, p[:, self.t].mean(axis=1)))
        newt = np.empty((self.t.shape[0], 2 * self.t.shape[1]),
                        dtype=self.t.dtype)
        newt[0, ::2] = self.t[0]
        newt[0, 1::2] = p.shape[1] + np.arange(self.t.shape[1])
        newt[1, ::2] = newt[0, 1::2]
        newt[1, 1::2] = self.t[1]
        # update fields
        self.p = newp
        self.facets = np.hstack(
            [self.facets,
             self.facets.shape[1] + np.arange(self.t.shape[1])[None, :]])
        self.t = newt
        self._build_mappings()

    def boundary_nodes(self):
        """Find the boundary nodes of the mesh."""
        _, counts = np.unique(self.t.flatten(), return_counts=True)
        return np.nonzero(counts == 1)[0]

    def interior_nodes(self):
        """Find the interior nodes of the mesh."""
        _, counts = np.unique(self.t.flatten(), return_counts=True)
        return np.nonzero(counts == 2)[0]

    def __mul__(self, other):
        """Tensor product mesh."""
        from .mesh2d.mesh_quad import MeshQuad
        return MeshQuad.init_tensor(self.p[0], other.p[0])

    def param(self):
        return np.max(np.abs(self.p[0, self.t[1, :]]
                             - self.p[0, self.t[0, :]]))

    def mapping(self):
        from skfem.mapping import MappingAffine
        return MappingAffine(self)

    @staticmethod
    def strip_extra_coordinates(p: ndarray) -> ndarray:
        return p[:, :1]

max_stars_repo_name: ahhuhtal/scikit-fem

max_issues_repo_name: ahhuhtal/scikit-fem

max_forks_repo_name: ahhuhtal/scikit-fem

reason: import numpy,from numpy


text:
import numpy as np

ads = np.zeros(4)
ads[0] = 4.031128e-09
ads[1] = 2.692970e-09
ads[2] = 1.594998e-09
ads[3] = 2.642999e-09

des = np.zeros(4)
des[0] = 6.592305e+11
des[1] = 1.844085e+11
des[2] = 1.161544e+10
des[3] = 2.293841e+09

dens = np.zeros(4)
#dens[0] = 1.817055e+19
#dens[1] = 3.604090e+18
#dens[2] = 1.820603e+18
#dens[3] = 8.678933e+17
dens[0] = 1.81694e+19
dens[1] = 3.60356e+18
dens[2] = 1.81807e+18
dens[3] = 8.57846e+17

t0 = 0.
T = 1e-12*10000

spec_init = np.zeros(4)
spec_init[0] = 0.;
spec_init[1] = 0.;
spec_init[2] = 0.;
spec_init[3] = 0.;

N = 100000
N_0 = 80*80

# calculate h (step-size)
h = (T - t0)/N

# define sol'n vectors
t = np.zeros(N+1)
spec1 = np.zeros(N+1,dtype=np.longdouble)
spec2 = np.zeros(N+1,dtype=np.longdouble)
spec3 = np.zeros(N+1,dtype=np.longdouble)
spec4 = np.zeros(N+1,dtype=np.longdouble)

t[0] = t0
spec1[0] = spec_init[0]
spec2[0] = spec_init[1]
spec3[0] = spec_init[2]
spec4[0] = spec_init[3]

f = lambda x,y,z,u: (1 - x - y - z - u)*ads[0]*dens[0] - des[0]*x
g = lambda x,y,z,u: (1 - x - y - z - u)*ads[1]*dens[1] - des[1]*y
p = lambda x,y,z,u: (1 - x - y - z - u)*ads[2]*dens[2] - des[2]*z
q = lambda x,y,z,u: (1 - x - y - z - u)*ads[3]*dens[3] - des[3]*u

F = [f,g,p,q]
num_F = len(F)

k = np.zeros([num_F,4],dtype=np.longdouble)

for i in range(0,N):
    for j in range(0,num_F):
        k[j,0] = h*F[j](spec1[i],spec2[i],spec3[i],spec4[i])
    for j in range(0,num_F):
        k[j,1] = h*F[j](spec1[i] + k[0,0]/2,spec2[i] + k[1,0]/2,spec3[i] + k[2,0]/2,spec4[i] + k[3,0]/2)
    for j in range(0,num_F):
        k[j,2] = h*F[j](spec1[i] + k[0,1]/2,spec2[i] + k[1,1]/2,spec3[i] + k[2,1]/2,spec4[i] + k[3,1]/2)
    for j in range(0,num_F):
        k[j,3] = h*F[j](spec1[i] + k[0,2]  ,spec2[i] + k[1,2]  ,spec3[i] + k[2,2]  ,spec4[i] + k[3,2])

    spec1[i+1] = spec1[i] + (1./6)*(k[0,0] + 2*k[0,1] + 2*k[0,2] + k[0,3])
    spec2[i+1] = spec2[i] + (1./6)*(k[1,0] + 2*k[1,1] + 2*k[1,2] + k[1,3])
    spec3[i+1] = spec3[i] + (1./6)*(k[2,0] + 2*k[2,1] + 2*k[2,2] + k[2,3])
    spec4[i+1] = spec4[i] + (1./6)*(k[3,0] + 2*k[3,1] + 2*k[3,2] + k[3,3])
    t[i+1] = t[i] + h

spec1 = N_0 * spec1
spec2 = N_0 * spec2
spec3 = N_0 * spec3
spec4 = N_0 * spec4

np.savetxt("res.rk4",np.transpose([t,spec1,spec2,spec3,spec4]))

max_stars_repo_name: PeculiarOvertones/FHDeX

max_issues_repo_name: PeculiarOvertones/FHDeX

max_forks_repo_name: PeculiarOvertones/FHDeX


text:
# -*- coding: utf-8 -*-
#
# Copyright 2019 Benjamin Kiessling
#
# 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.
"""
kraken.blla
~~~~~~~~~~~~~~

Trainable baseline layout analysis tools for kraken
"""

import torch
import logging
import numpy as np
import pkg_resources
import torch.nn.functional as F
import torchvision.transforms as tf

from typing import Optional, Dict, Callable

from kraken.lib import vgsl, dataset
from kraken.lib.util import pil2array, is_bitonal, get_im_str
from kraken.lib.exceptions import KrakenInputException
from kraken.lib.segmentation import (polygonal_reading_order,
                                     vectorize_lines, vectorize_regions,
                                     scale_polygonal_lines,
                                     calculate_polygonal_environment,
                                     scale_regions)

__all__ = ['segment']

logger = logging.getLogger(__name__)

def segment(im,
            text_direction: str = 'horizontal-lr',
            mask: Optional[np.array] = None,
            reading_order_fn: Callable = polygonal_reading_order,
            model=None,
            device: str = 'cpu'):
    """
    Segments a page into text lines using the baseline segmenter.

    Segments a page into text lines and returns the polyline formed by each
    baseline and their estimated environment.

    Args:
        im (PIL.Image): An RGB image.
        text_direction (str): Ignored by the segmenter but kept for
                              serialization.
        mask (PIL.Image): A bi-level mask image of the same size as `im` where
                          0-valued regions are ignored for segmentation
                          purposes. Disables column detection.
        reading_order_fn (function): Function to determine the reading order.
                                     Has to accept a list of tuples (baselines,
                                     polygon) and a text direction (`lr` or
                                     `rl`).
        model (vgsl.TorchVGSLModel): A TorchVGSLModel containing a segmentation
                                     model. If none is given a default model
                                     will be loaded.
        device (str or torch.Device): The target device to run the neural
                                      network on.

    Returns:
        {'text_direction': '$dir',
         'type': 'baseline',
         'lines': [
            {'baseline': [[x0, y0], [x1, y1], ..., [x_n, y_n]], 'boundary': [[x0, y0, x1, y1], ... [x_m, y_m]]},
            {'baseline': [[x0, ...]], 'boundary': [[x0, ...]]}
          ]
          'regions': [
            {'region': [[x0, y0], [x1, y1], ..., [x_n, y_n]], 'type': 'image'},
            {'region': [[x0, ...]], 'type': 'text'}
          ]
        }: A dictionary containing the text direction and under the key 'lines'
        a list of reading order sorted baselines (polylines) and their
        respective polygonal boundaries. The last and first point of each
        boundary polygon is connected.

    Raises:
        KrakenInputException if the input image is not binarized or the text
        direction is invalid.
    """
    im_str = get_im_str(im)
    logger.info(f'Segmenting {im_str}')

    if model is None:
        logger.info('No segmentation model given. Loading default model.')
        model = vgsl.TorchVGSLModel.load_model(pkg_resources.resource_filename(__name__, 'blla.mlmodel'))

    if model.one_channel_mode == '1' and not is_bitonal(im):
        logger.warning('Running binary model on non-binary input image '
                       '(mode {}). This will result in severely degraded '
                       'performance'.format(im.mode))

    model.eval()
    model.to(device)

    if mask:
        if mask.mode != '1' and not is_bitonal(mask):
            logger.error('Mask is not bitonal')
            raise KrakenInputException('Mask is not bitonal')
        mask = mask.convert('1')
        if mask.size != im.size:
            logger.error('Mask size {mask.size} doesn\'t match image size {im.size}')
            raise KrakenInputException('Mask size {mask.size} doesn\'t match image size {im.size}')
        logger.info('Masking enabled in segmenter.')
        mask = pil2array(mask)

    batch, channels, height, width = model.input
    transforms = dataset.generate_input_transforms(batch, height, width, channels, 0, valid_norm=False)
    res_tf = tf.Compose(transforms.transforms[:3])
    scal_im = res_tf(im).convert('L')

    with torch.no_grad():
        logger.debug('Running network forward pass')
        o, _ = model.nn(transforms(im).unsqueeze(0).to(device))
    logger.debug('Upsampling network output')
    o = F.interpolate(o, size=scal_im.size[::-1])
    o = o.squeeze().cpu().numpy()
    scale = np.divide(im.size, o.shape[:0:-1])
    # postprocessing
    cls_map = model.user_metadata['class_mapping']
    st_sep = cls_map['aux']['_start_separator']
    end_sep = cls_map['aux']['_end_separator']

    logger.info('Vectorizing baselines')
    baselines = []
    regions = {}
    for bl_type, idx in cls_map['baselines'].items():
        logger.debug(f'Vectorizing lines of type {bl_type}')
        baselines.extend([(bl_type,x) for x in vectorize_lines(o[(st_sep, end_sep, idx), :, :])])
    logger.info('Vectorizing regions')
    for region_type, idx in cls_map['regions'].items():
        logger.debug(f'Vectorizing lines of type {bl_type}')
        regions[region_type] = vectorize_regions(o[idx])
    logger.debug('Polygonizing lines')
    lines = list(filter(lambda x: x[2] is not None, zip([x[0] for x in baselines],
                                                        [x[1] for x in baselines],
                                                        calculate_polygonal_environment(scal_im, [x[1] for x in baselines]))))
    logger.debug('Scaling vectorized lines')
    sc = scale_polygonal_lines([x[1:] for x in lines], scale)
    lines = list(zip([x[0] for x in lines], [x[0] for x in sc], [x[1] for x in sc]))
    logger.debug('Scaling vectorized regions')
    for reg_id, regs in regions.items():
        regions[reg_id] = scale_regions(regs, scale)
    logger.debug('Reordering baselines')
    order_regs = []
    for regs in regions.values():
        order_regs.extend(regs)
    lines = reading_order_fn(lines=lines, regions=order_regs, text_direction=text_direction[-2:])

    if 'class_mapping' in model.user_metadata and len(model.user_metadata['class_mapping']['baselines']) > 1:
        script_detection = True
    else:
        script_detection = False

    return {'text_direction': text_direction,
            'type': 'baselines',
            'lines': [{'script': bl_type, 'baseline': bl, 'boundary': pl} for bl_type, bl, pl in lines],
            'regions': regions,
            'script_detection': script_detection}

max_stars_repo_name: eighttails/kraken

max_issues_repo_name: eighttails/kraken

max_forks_repo_name: eighttails/kraken


text:
#!/usr/bin/env python3.6

import os
import numpy as np
import sys
try:
  import torch
except ImportError:
    pass
from easypbr  import *
from dataloaders import *
# np.set_printoptions(threshold=sys.maxsize)
np.random.seed(0)

config_file="test_loader.cfg"

config_path=os.path.join( os.path.dirname( os.path.realpath(__file__) ) , '../config', config_file)
view=Viewer.create(config_path) #first because it needs to init context

def test_volref():
    loader=DataLoaderVolRef(config_path)
    loader.start()

    while True:
        if(loader.has_data()  ):

            first=True

            #volref
            # print("got frame")
            frame_color=loader.get_color_frame()
            frame_depth=loader.get_depth_frame()

            Gui.show(frame_color.rgb_32f, "rgb")

            rgb_with_valid_depth=frame_color.rgb_with_valid_depth(frame_depth)
            Gui.show(rgb_with_valid_depth, "rgb_valid")

            frustum_mesh=frame_depth.create_frustum_mesh(0.02)
            frustum_mesh.m_vis.m_line_width=3
            frustum_name="frustum"
            Scene.show(frustum_mesh, frustum_name)
            # Scene.show(frustum_mesh, "frustum_"+str(frame_color.frame_idx) )

            cloud=frame_depth.depth2world_xyz_mesh()
            frame_color.assign_color(cloud) #project the cloud into this frame and creates a color matrix for it
            Scene.show(cloud, "cloud")

        if loader.is_finished():
            # print("resetting")
            loader.reset()

        view.update()

def test_img():
    loader=DataLoaderImg(config_path)
    loader.start()

    while True:
        for cam_idx in range(loader.get_nr_cams()):
            if(loader.has_data_for_cam(cam_idx)):
                # print("data")
                frame=loader.get_frame_for_cam(cam_idx)
                Gui.show(frame.rgb_32f, "rgb")
                #get tensor
                rgb_tensor=frame.rgb2tensor()
                rgb_tensor=rgb_tensor.to("cuda")

        view.update()

def test_semantickitti():
    loader=DataLoaderSemanticKitti(config_path)
    loader.start()

    while True:
        if(loader.has_data()):
            cloud=loader.get_cloud()
            Scene.show(cloud, "cloud")

            # if cloud.V.size(0)==125620:
                # print("found")
        if loader.is_finished():
            loader.reset()

        view.update()

def test_scannet():
    loader=DataLoaderScanNet(config_path)
    loader.start()

    while True:
        if(loader.has_data()):
            cloud=loader.get_cloud()
            Scene.show(cloud, "cloud")

        view.update()

def test_stanford3dscene():
    loader=DataLoaderStanford3DScene(config_path)
    loader.start()

    while True:
        if(loader.has_data() ):

            print("got frame")
            frame_color=loader.get_color_frame()
            frame_depth=loader.get_depth_frame()

            Gui.show(frame_color.rgb_32f, "rgb")

            rgb_with_valid_depth=frame_color.rgb_with_valid_depth(frame_depth)
            Gui.show(rgb_with_valid_depth, "rgb_valid")

            frustum_mesh=frame_color.create_frustum_mesh(0.1)
            frustum_mesh.m_vis.m_line_width=3
            frustum_name="frustum"
            Scene.show(frustum_mesh, frustum_name)

            cloud=frame_depth.backproject_depth()
            frame_color.assign_color(cloud) #project the cloud into this frame and creates a color matrix for it
            # Scene.show(cloud, "cloud"+str(nr_cloud))
            Scene.show(cloud, "cloud")

            #show look dir
            look_dir_mesh=Mesh()
            look_dir=frame_color.look_dir()
            look_dir_mesh.V=[
                look_dir
            ]
            look_dir_mesh.m_vis.m_show_points=True
            Scene.show(look_dir_mesh, "look_dir_mesh")

        view.update()

def test_shapenet_img():
    loader=DataLoaderShapeNetImg(config_path)

    i=0

    while True:
        if(loader.finished_reading_scene() ):
            frame=loader.get_random_frame()

            if i%1==0:
                loader.start_reading_next_scene()

            if frame.is_shell:
                frame.load_images()

            Gui.show(frame.rgb_32f, "rgb")
            Gui.show(frame.mask, "mask")
            Gui.show(frame.depth, "depth")
            frustum=frame.create_frustum_mesh(0.1)
            # Scene.show(frustum, "frustum"+ str(frame.frame_idx) )
            Scene.show(frustum, "frustum" )

            # cloud=frame.depth2world_xyz_mesh()
            # cloud=frame.assign_color(cloud)
            # Scene.show(cloud, "cloud")

            i+=1

        if loader.is_finished():
            print("resetting")
            loader.reset()

        view.update()

def test_nerf():
    loader=DataLoaderNerf(config_path)
    loader.set_mode_train()
    # loader.set_mode_test()
    loader.start()

    while True:
        if(loader.has_data() ):

            # print("got frame")
            frame=loader.get_next_frame()
            # print("frame width and height is ", frame.width, " ", frame.height)

            Gui.show(frame.rgb_32f, "rgb")

            frustum_mesh=frame.create_frustum_mesh(0.02)
            frustum_mesh.m_vis.m_line_width=1
            Scene.show(frustum_mesh, "frustum_"+str(frame.frame_idx) )

            # cloud=frame.depth2world_xyz_mesh()
            # Scene.show(cloud, "cloud")

        if loader.is_finished():
            print("resetting")
            loader.reset()
            print("scene scale is ", Scene.get_scale())

        view.update()

def test_phenorob():
    loader=DataLoaderPhenorob(config_path)
    loader.start()
    loader.set_do_augmentation(True)

    while True:
        if(loader.has_data() ):
            # print("has data")

            cloud=loader.get_cloud()

            Scene.show(cloud, "cloud" )

        if loader.is_finished():
            # print("resetting")
            loader.reset()

        view.update()

def test_cloud_ros():
    bag=RosBagPlayer.create(config_path)
    loader=DataLoaderCloudRos(config_path)

    while loader.is_loader_thread_alive():
        if(loader.has_data() ):
            cloud=loader.get_cloud()

            Scene.show(cloud, "cloud" )

        view.update()

def test_colmap():
    loader=DataLoaderColmap(config_path)
    loader.start()

    while True:
        if(loader.has_data() ):

            # print("got frame")
            frame=loader.get_next_frame()
            # print("loaded frame with width and height ", frame.width, " ", frame.height)

            Gui.show(frame.rgb_32f, "rgb")

            # print("frame k is ", frame.K)
            frustum_mesh=frame.create_frustum_mesh(0.02)
            frustum_mesh.m_vis.m_line_width=1
            Scene.show(frustum_mesh, "frustum_"+str(frame.frame_idx) )

            # cloud=frame.depth2world_xyz_mesh()
            # Scene.show(cloud, "cloud")

        if loader.is_finished():
            print("resetting")
            loader.reset()
            print("scene scale is ", Scene.get_scale())

        view.update()

def test_easypbr():
    loader=DataLoaderEasyPBR(config_path)
    loader.set_mode_train()
    # loader.set_mode_test()
    loader.start()

    while True:
        if(loader.has_data() ):

            # print("got frame")
            frame=loader.get_next_frame()
            # print("frame width and height is ", frame.width, " ", frame.height)

            # if frame.frame_idx==20:
            Gui.show(frame.rgb_32f, "rgb")

            # frustum_mesh=frame.create_frustum_mesh(0.02)
            frustum_mesh=frame.create_frustum_mesh(0.2)
            frustum_mesh.m_vis.m_line_width=1
            if frame.frame_idx==20:
                frustum_mesh.m_vis.m_line_color=[0.0, 1.0, 0.0]
            Scene.show(frustum_mesh, "frustum_"+str(frame.frame_idx) )
            # mesh=Scene.get_mesh_with_name("frustum_"+str(frame.frame_idx))
            # frustum_mesh.m_vis.m_line_width=1
            # frustum_mesh.m_vis.m_line_color=[0.0, 1.0, 0.0]

            # #make all frustums red except the current one
            # for i in range(200):
            #     if i==frame.frame_idx:
            #         continue
            #     if Scene.does_mesh_with_name_exist("frustum_"+str(i)):
            #         mesh=Scene.get_mesh_with_name("frustum_"+str(i))
            #         frustum_mesh.m_vis.m_line_width=1
            #         frustum_mesh.m_vis.m_line_color=[1.0, 0.0, 0.0]

            # cloud=frame.depth2world_xyz_mesh()
            # Scene.show(cloud, "cloud")

        if loader.is_finished():
            print("resetting")
            loader.reset()
            print("scene scale is ", Scene.get_scale())

        view.update()

def test_srn():
    loader=DataLoaderSRN(config_path)
    loader.set_mode_train()
    # loader.set_mode_test()
    loader.start()

    i=0

    while True:
        if(loader.finished_reading_scene() ):
            frame=loader.get_random_frame()

            if i%1==0:
                loader.start_reading_next_scene()

            if frame.is_shell:
                frame.load_images()

            Gui.show(frame.rgb_32f, "rgb")
            frustum=frame.create_frustum_mesh(0.02)
            Scene.show(frustum, "frustum"+ str(frame.frame_idx) )
            # Scene.show(frustum, "frustum" )

            # cloud=frame.depth2world_xyz_mesh()
            # cloud=frame.assign_color(cloud)
            # Scene.show(cloud, "cloud")

            i+=1

        if loader.is_finished():
            print("resetting")
            loader.reset()

        view.update()

def test_dtu():
    loader=DataLoaderDTU(config_path)
    loader.set_mode_train()
    # loader.set_mode_validation() #test set actually doesnt exist and we actually use the validation one
    loader.start()

    i=0

    while True:
        if(loader.finished_reading_scene() ):
            frame=loader.get_random_frame()

            # if i%1==0:
                # loader.start_reading_next_scene()

            print("frame rgb path is ", frame.rgb_path)

            if frame.is_shell:
                frame.load_images()

            Gui.show(frame.rgb_32f, "rgb")
            frustum=frame.create_frustum_mesh(0.02)
            Scene.show(frustum, "frustum"+ str(frame.frame_idx) )
            # Scene.show(frustum, "frustum" )

            loader.start_reading_next_scene()
            while True:
                if(loader.finished_reading_scene()):
                    break

            i+=1

        if loader.is_finished():
            print("resetting")
            loader.reset()

        view.update()

def test_deep_voxels():
    loader=DataLoaderDeepVoxels(config_path)
    loader.set_mode_train()
    # loader.set_mode_test()
    loader.start()

    while True:
        if(loader.has_data() ):

            # print("got frame")
            frame=loader.get_next_frame()
            # print("frame width and height is ", frame.width, " ", frame.height)

            Gui.show(frame.rgb_32f, "rgb")

            frustum_mesh=frame.create_frustum_mesh(0.02)
            frustum_mesh.m_vis.m_line_width=1
            Scene.show(frustum_mesh, "frustum_"+str(frame.frame_idx) )

            # cloud=frame.depth2world_xyz_mesh()
            # Scene.show(cloud, "cloud")

        if loader.is_finished():
            print("resetting")
            loader.reset()
            print("scene scale is ", Scene.get_scale())

        view.update()

def test_llff():
    loader=DataLoaderLLFF(config_path)
    loader.set_mode_train()
    loader.start()

    while True:
        if(loader.has_data() ):

            # print("got frame")
            frame=loader.get_next_frame()
            # print("loaded frame with width and height ", frame.width, " ", frame.height)

            Gui.show(frame.rgb_32f, "rgb")

            # print("frame k is ", frame.K)
            frustum_mesh=frame.create_frustum_mesh(0.02)
            frustum_mesh.m_vis.m_line_width=1
            Scene.show(frustum_mesh, "frustum_"+str(frame.frame_idx) )

            # cloud=frame.depth2world_xyz_mesh()
            # Scene.show(cloud, "cloud")

            near = frame.get_extra_field_float("min_near")
            far = frame.get_extra_field_float("max_far")
            print("near is ", near)
            print("far is ", far)
            # frame.add_extra_field("test", 54.0)
            # roundback= frame.get_extra_field_float("test")
            # print("roundback is ", roundback)

        if loader.is_finished():
            print("resetting")
            loader.reset()
            print("scene scale is ", Scene.get_scale())

        view.update()

def test_blender_fb():
    loader=DataLoaderBlenderFB(config_path)
    # loader.set_mode_train()
    # loader.set_mode_test()
    loader.start()

    while True:
        if(loader.has_data() ):

            # print("got frame")
            frame=loader.get_next_frame()
            # print("frame width and height is ", frame.width, " ", frame.height)

            Gui.show(frame.rgb_32f, "rgb")
            # Gui.show(frame.rgb_8u, "rgb8u")

            frustum_mesh=frame.create_frustum_mesh(0.05)
            frustum_mesh.m_vis.m_line_width=1
            # frustum_mesh.m_is_dirty=True
            Scene.show(frustum_mesh, "frustum_"+str(frame.frame_idx) )

            # cloud=frame.depth2world_xyz_mesh()
            # Scene.show(cloud, "cloud")

            if frame.has_extra_field("orientation_mat"):
                orientation_mat = frame.get_extra_field_mat("orientation_mat")
                variance_mat = frame.get_extra_field_mat("variance_mat")
                Gui.show(orientation_mat,"orientation_mat")
                Gui.show(variance_mat,"variance_mat")

        if loader.is_finished():
            # print("resetting")
            loader.reset()
            # print("scene scale is ", Scene.get_scale())

        view.update()

def test_usc_hair():
    loader=DataLoaderUSCHair(config_path)
    loader.set_mode_train()
    # loader.set_mode_test()
    loader.start()

    while True:
        if(loader.has_data()):
            # print("HAS DATA--------")
            # pass
            # hair=loader.get_hair()
            # cloud=hair.full_hair_cloud
            cloud=loader.get_hair().full_hair_cloud.clone()
            Scene.show(cloud, "cloud")

            #show also the head
            # head=loader.get_mesh_head()
            # Scene.show(head, "head")

            #show also the scalp
            scalp=loader.get_mesh_scalp()
            Scene.show(scalp, "scalp")

            if cloud.has_extra_field("strand_idx"):
                strand_idx=cloud.get_extra_field_matrixXi("strand_idx")
                # print("strand idx", strand_idx)

            if cloud.has_extra_field("uv_roots"):
                uv_roots=cloud.get_extra_field_matrixXd("uv_roots")
                # print("uv_roots", uv_roots)

        if loader.is_finished():
            loader.reset()

        view.update()

# test_volref()
# test_img()
# test_img_ros()
# test_cloud_ros()
# test_semantickitti()
# test_scannet()
# test_stanford3dscene()
# test_shapenet_img()
# test_nerf()
# test_phenorob()
# test_colmap()
# test_easypbr()
# test_srn()
# test_dtu()
# test_deep_voxels()
# test_llff()
# test_blender_fb()
test_usc_hair()


text:
import os
import pandas as pd
from plots import *
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np

######################
# DONE: LOAD DATASET #
######################

project_path = os.sep.join(os.getcwd().split(os.sep)[1:-1])

data_root_path = os.path.join(os.sep, project_path, 'sprint_03_data_analysis', 'output')

file_csv_output_01 = 'output_01_data_01_AM_mun_especie.csv'

path_sprint_03_data_analysis_output_01_csv = os.path.join(data_root_path, file_csv_output_01)
df = pd.read_csv(path_sprint_03_data_analysis_output_01_csv, delimiter=';')

del df['Unnamed: 0']

df = df[df['Estado'] != 'Total']

print(df.info())

"""
RangeIndex: 9101 entries, 0 to 9100
Data columns (total 10 columns):
 #   Column      Non-Null Count  Dtype
---  ------      --------------  -----
 0   Unnamed: 0  9101 non-null   int64
 1   CD          9101 non-null   object
 2   Municipio   9101 non-null   object
 3   Falciparum  9101 non-null   float64
 4   Mista       9101 non-null   float64
 5   Vivax       9101 non-null   float64
 6   Malarie     9101 non-null   int64
 7   Ovale       9101 non-null   int64
 8   Ano         9101 non-null   int64
 9   Estado      9101 non-null   object
dtypes: float64(3), int64(4), object(3)
memory usage: 711.1+ KB
"""
#############################
# TODO: DATA VISUALIZATIONS #
#############################

##############################
# Data analysis exploitation #
##############################

df.describe()

print(df['Vivax'].sum()) # 795600.3940000001

print(df['Falciparum'].sum()) # 402047.081

print(df['Malarie'].sum()) # 3674

int(df['Ovale'].sum()) # 102

"""
Vivax > Falciparum > Malarie > Ovale
"""

####################################
# DONE: Distribution plots         #
# IMPOSSIBLE WITH CATEGORICAL DATA #
####################################

path_output = os.path.join(os.getcwd(), 'output')

distplot(df, 'Falciparum', kde=True).savefig(os.path.join(path_output, 'output_01.png'))

distplot(df, 'Vivax', kde=True).savefig(os.path.join(path_output, 'output_02.png'))

distplot(df, 'Malarie', kde=True).savefig(os.path.join(path_output, 'output_03.png'))

distplot(df, 'Ovale', kde=True).savefig(os.path.join(path_output, 'output_04.png'))

jointplot(df, 'Estado', 'Falciparum').savefig(os.path.join(path_output, 'output_05.png'))

jointplot(df, 'Estado', 'Vivax').savefig(os.path.join(path_output, 'output_06.png'))

jointplot(df, 'Estado', 'Malarie').savefig(os.path.join(path_output, 'output_07.png'))

jointplot(df, 'Estado', 'Ovale').savefig(os.path.join(path_output, 'output_08.png'))

jointplot(df, 'Ano', 'Falciparum').savefig(os.path.join(path_output, 'output_09.png'))

jointplot(df, 'Ano', 'Vivax').savefig(os.path.join(path_output, 'output_10.png'))

jointplot(df, 'Ano', 'Malarie').savefig(os.path.join(path_output, 'output_11.png'))

jointplot(df, 'Ano', 'Ovale').savefig(os.path.join(path_output, 'output_12.png'))

pairplot(df).savefig(os.path.join(path_output, 'output_13.png'))

sns.set(rc={'figure.figsize': (10, 10)})

pairplot(df[['Falciparum', 'Vivax', 'Estado']], hue='Estado').savefig(os.path.join(path_output, 'output_14.png'))

pairplot(df[['Falciparum', 'Vivax', 'Estado', 'Ano']], hue='Estado').savefig(os.path.join(path_output, 'output_15.png'))

pairplot(df[['Falciparum', 'Vivax', 'Estado', 'Ano']], hue='Estado', kind='kde').savefig(os.path.join(path_output, 'output_16.png'))

pairplot(df[['Falciparum', 'Vivax', 'Estado']], hue='Estado', kind='kde').savefig(os.path.join(path_output, 'output_17.png'))

pairplot(df[['Falciparum', 'Vivax', 'Estado']], hue='Estado', kind='reg').savefig(os.path.join(path_output, 'output_18.png'))

pairplot(df[['Falciparum', 'Vivax', 'Estado']], hue='Estado', kind='reg').savefig(os.path.join(path_output, 'output_19.png'))

pairplot(df[['Falciparum', 'Vivax', 'Estado', 'Ano']][df['Ano'] >= 2020], hue='Estado', kind='reg').savefig(os.path.join(path_output, 'output_20.png'))

pairplot(df[['Falciparum', 'Vivax', 'Estado', 'Ano']][df['Ano'] >= 2020], hue='Estado', kind='hex').savefig(os.path.join(path_output, 'output_21.png'))

pairplot(df[['Falciparum', 'Vivax', 'Ano']][df['Ano'] >= 2020], hue='Ano', kind='scatter', diag_kind='kde').savefig(os.path.join(path_output, 'output_22.png'))

pairplot(df[['Falciparum', 'Vivax', 'Ano']][df['Ano'] >= 2020], hue='Ano', kind='scatter', diag_kind='hist').savefig(os.path.join(path_output, 'output_23.png'))

pairplot(df[['Falciparum', 'Vivax', 'Ano']][df['Ano'] >= 2020], hue='Ano', kind='kde', diag_kind='kde').savefig(os.path.join(path_output, 'output_24.png'))

pairplot(df[['Falciparum', 'Vivax', 'Ano']][df['Ano'] >= 2020], hue='Ano', kind='reg', diag_kind='kde').savefig(os.path.join(path_output, 'output_25.png'))

##########################
# DONE: CATEGORICAL PLOT #
##########################

barplot(df, x='Estado', y='Falciparum').figure.savefig(os.path.join(path_output, 'output_26.png'))

barplot(df, x='Estado', y='Vivax').figure.savefig(os.path.join(path_output, 'output_27.png'))

barplot(df, x='Ano', y='Falciparum').figure.savefig(os.path.join(path_output, 'output_28.png'))

barplot(df, x='Ano', y='Vivax').figure.savefig(os.path.join(path_output, 'output_29.png'))

sns.set(rc={'figure.figsize': (30, 10)})

barplot(df[df['Ano'] >= 2010], x='Estado', y='Falciparum', hue='Ano').figure.savefig(os.path.join(path_output, 'output_30.png'))

barplot(df[df['Ano'] >= 2010], x='Estado', y='Vivax', hue='Ano').figure.savefig(os.path.join(path_output, 'output_31.png'))

barplot(df[df['Ano'] >= 2010], x='Ano', y='Falciparum', hue='Estado').figure.savefig(os.path.join(path_output, 'output_32.png'))

barplot(df[df['Ano'] >= 2010], x='Ano', y='Vivax', hue='Estado').figure.savefig(os.path.join(path_output, 'output_33.png'))

sns.set(rc={'figure.figsize': (10, 10)})

countplot(df, x='Estado').set_title('Municípios por estado').figure.savefig(os.path.join(path_output, 'output_34.png'))

boxplot(df, x='Estado', y='Falciparum').set_title('Falciparum por estado').figure.savefig(os.path.join(path_output, 'output_35.png'))

boxplot(df, x='Ano', y='Falciparum').set_title('Falciparum por ano').figure.savefig(os.path.join(path_output, 'output_36.png'))

boxplot(df, x='Estado', y='Vivax').set_title('Vivax por estado').figure.savefig(os.path.join(path_output, 'output_37.png'))

boxplot(df, x='Ano', y='Vivax').set_title('Vivax por ano').figure.savefig(os.path.join(path_output, 'output_38.png'))

sns.set(rc={'figure.figsize': (30, 10)})

boxplot(df, x='Estado', y='Falciparum', hue='Ano').set_title('Falciparum por Estado e ano').figure.savefig(os.path.join(path_output, 'output_39.png'))

boxplot(df, x='Ano', y='Falciparum', hue='Estado').set_title('Falciparum por ano e Estado').figure.savefig(os.path.join(path_output, 'output_40.png'))

boxplot(df, x='Estado', y='Vivax', hue='Ano').set_title('Vivax por Estado e ano').figure.savefig(os.path.join(path_output, 'output_41.png'))

boxplot(df, x='Ano', y='Vivax', hue='Estado').set_title('Vivax por ano e Estado').figure.savefig(os.path.join(path_output, 'output_42.png'))

sns.set(rc={'figure.figsize': (20, 10)})

violinplot(df, x='Estado', y='Falciparum').set_title('Falciparum por estado').figure.savefig(os.path.join(path_output, 'output_43.png'))

violinplot(df, x='Ano', y='Falciparum').set_title('Falciparum por ano').figure.savefig(os.path.join(path_output, 'output_44.png'))

violinplot(df, x='Estado', y='Vivax').set_title('Vivax por estado').figure.savefig(os.path.join(path_output, 'output_45.png'))

violinplot(df, x='Ano', y='Vivax').set_title('Vivax por ano').figure.savefig(os.path.join(path_output, 'output_46.png'))

sns.set(rc={'figure.figsize': (30, 10)})

stripplot(df, x='Estado', y='Falciparum', hue='Ano', dodge=True).set_title('Falciparum por Estado e ano').figure.savefig(os.path.join(path_output, 'output_47.png'))

stripplot(df, x='Ano', y='Falciparum', hue='Estado', dodge=True).set_title('Falciparum por ano e Estado').figure.savefig(os.path.join(path_output, 'output_48.png'))

stripplot(df, x='Estado', y='Vivax', hue='Ano', dodge=True).set_title('Vivax por Estado e ano').figure.savefig(os.path.join(path_output, 'output_49.png'))

stripplot(df, x='Ano', y='Vivax', hue='Estado', dodge=True).set_title('Vivax por ano e Estado').figure.savefig(os.path.join(path_output, 'output_50.png'))

catplot(df, x='Ano', y='Falciparum', col='Estado', kind='bar', aspect=1.7, col_wrap=2).savefig(os.path.join(path_output, 'output_51.png'))

catplot(df, x='Ano', y='Falciparum', col='Estado', kind='box', aspect=1.7, col_wrap=2).savefig(os.path.join(path_output, 'output_52.png'))

catplot(df, x='Ano', y='Falciparum', col='Estado', kind='strip', aspect=1.7, col_wrap=2).savefig(os.path.join(path_output, 'output_53.png'))

catplot(df, x='Ano', y='Vivax', col='Estado', kind='bar', aspect=1.7, col_wrap=2).savefig(os.path.join(path_output, 'output_54.png'))

catplot(df, x='Ano', y='Vivax', col='Estado', kind='box', aspect=1.7, col_wrap=2).savefig(os.path.join(path_output, 'output_55.png'))

catplot(df, x='Ano', y='Vivax', col='Estado', kind='strip', aspect=1.7, col_wrap=2).savefig(os.path.join(path_output, 'output_56.png'))

#########################
# DONE: MATRICIAL PLOTS #
#########################

sns.set(rc={'figure.figsize': (20, 10)})

heatmap(df.pivot_table(values='Falciparum', index='Estado', columns='Ano'), cmap='magma').set_title('Falciparum').figure.savefig(os.path.join(path_output, 'output_57.png'))

heatmap(df.pivot_table(values='Vivax', index='Estado', columns='Ano'), cmap='magma').set_title('Vivax').figure.savefig(os.path.join(path_output, 'output_58.png'))

clustermap(df.pivot_table(values='Falciparum', index='Estado', columns='Ano'), cmap='magma').tsavefig(os.path.join(path_output, 'output_59.png'))

clustermap(df.pivot_table(values='Vivax', index='Estado', columns='Ano'), cmap='magma').savefig(os.path.join(path_output, 'output_60.png'))

####################
# REGRESSION PLOTS #
####################

sns.lmplot(data=df, x='Ano', y='Falciparum', col='Estado', aspect=1.7, col_wrap=2, palette='Set1').savefig(os.path.join(path_output, 'output_61.png'))

sns.lmplot(data=df, x='Ano', y='Vivax', col='Estado', aspect=1.7, col_wrap=2, palette='Set1').savefig(os.path.join(path_output, 'output_62.png'))

lmplot(data=df[df['Ano'] >= 2020], x='Ano', y='Falciparum', col='Estado', scatter=False, height=3, aspect=1, col_wrap=3, palette='Set1').set(xticks=[2020, 2021], yticks=[0, 100, 200, 300]).savefig(os.path.join(path_output, 'output_63.png'))

lmplot(data=df[df['Ano'] >= 2020], x='Ano', y='Vivax', col='Estado', scatter=False, height=3, aspect=1, col_wrap=3, palette='Set1').set(xticks=[2020, 2021], yticks=[0, 100, 200, 300]).savefig(os.path.join(path_output, 'output_64.png'))


text:
#!/usr/bin/env python3

from dataclasses import dataclass, field
import array
import numpy as np
import pdb
import copy

@dataclass
class PythonMsg:
    '''
    Base class for python messages. Intention is that fields cannot be changed accidentally,
    e.g. if you try "state.xk = 10", it will throw an error if no field "xk" exists.
    This helps avoid typos. If you really need to add a field use "object.__setattr__(state,'xk',10)"
    Dataclasses autogenerate a constructor, e.g. a dataclass fields x,y,z can be instantiated
    "pos = Position(x = 1, y = 2, z = 3)" without you having to write the __init__() function, just add the decorator
    Together it is hoped that these provide useful tools and safety measures for storing and moving data around by name,
    rather than by magic indices in an array, e.g. q = [9, 4.5, 8829] vs. q.x = 10, q.y = 9, q.z = 16
    '''
    def __setattr__(self, key, value):
        '''
        Overloads default attribute-setting functionality to avoid creating new fields that don't already exist
        This exists to avoid hard-to-debug errors from accidentally adding new fields instead of modifying existing ones

        To avoid this, use:
        object.__setattr__(instance, key, value)
        ONLY when absolutely necessary.
        '''
        if not hasattr(self, key):
            raise TypeError('Cannot add new field "%s" to frozen class %s' % (key, self))
        else:
            object.__setattr__(self, key, value)

    def print(self, depth=0, name=None):
        '''
        default __str__ method is not easy to read, especially for nested classes.
        This is easier to read but much longer

        Will not work with "from_str" method.
        '''
        print_str = ''
        for j in range(depth): print_str += '  '
        if name:
            print_str += name + ' (' + type(self).__name__ + '):\n'
        else:
            print_str += type(self).__name__ + ':\n'
        for key in vars(self):
            val = self.__getattribute__(key)
            if isinstance(val, PythonMsg):
                print_str += val.print(depth=depth + 1, name=key)
            else:
                for j in range(depth + 1): print_str += '  '
                print_str += str(key) + '=' + str(val)
                print_str += '\n'

        if depth == 0:
            print(print_str)
        else:
            return print_str

    def from_str(self, string_rep):
        '''
        inverts dataclass.__str__() method generated for this class so you can unpack objects sent via text (e.g. through multiprocessing.Queue)
        '''
        val_str_index = 0
        for key in vars(self):
            val_str_index = string_rep.find(key + '=', val_str_index) + len(key) + 1  # add 1 for the '=' sign
            value_substr = string_rep[val_str_index: string_rep.find(',',
                                                                     val_str_index)]  # (thomasfork) - this should work as long as there are no string entries with commas

            if '\'' in value_substr:  # strings are put in quotes
                self.__setattr__(key, value_substr[1:-1])
            if 'None' in value_substr:
                self.__setattr__(key, None)
            else:
                self.__setattr__(key, float(value_substr))

    def copy(self):
        return copy.deepcopy(self)

@dataclass
class Position(PythonMsg):
    x: float = field(default=0)
    y: float = field(default=0)
    z: float = field(default=0)

    # def xdot(self, q: 'OrientationQuaternion', v: 'BodyLinearVelocity') -> 'Position':
    #     xdot = Position()
    #     xdot.x = (1 - 2 * q.qj ** 2 - 2 * q.qk ** 2) * v.v1 + 2 * (q.qi * q.qj - q.qk * q.qr) * v.v2 + 2 * (
    #                 q.qi * q.qk + q.qj * q.qr) * v.v3
    #     xdot.y = (1 - 2 * q.qk ** 2 - 2 * q.qi ** 2) * v.v2 + 2 * (q.qj * q.qk - q.qi * q.qr) * v.v3 + 2 * (
    #                 q.qj * q.qi + q.qk * q.qr) * v.v1
    #     xdot.z = (1 - 2 * q.qi ** 2 - 2 * q.qj ** 2) * v.v3 + 2 * (q.qk * q.qi - q.qj * q.qr) * v.v1 + 2 * (
    #                 q.qk * q.qj + q.qi * q.qr) * v.v2
    #     return xdot

@dataclass
class VehicleActuation(PythonMsg):
    t: float = field(default = 0)

    u_a: float = field(default = 0)
    u_steer: float = field(default = 0)

    def __str__(self):
        return 't:{self.t}, u_a:{self.u_a}, u_steer:{self.u_steer}'.format(self=self)

@dataclass
class TrackLookahead(PythonMsg):
    '''
    Local track information ahead of the vehicle (curvature)
    '''
    t: float = field(default = None)    # time in seconds

    l: float = field(default = None)    # length of lookahead in meters
    dl: float = field(default = None)   # discretization step-size of the lookahead
    n: int = field(default = None)    # length of lookahead in array entries

    # TODO Add field for segmented lookahead?
    curvature: array.array = field(default = None)  # the curvature lookahead

    def __post_init__(self):
        if self.l is None: self.l = 1
        if self.dl is None: self.dl = 0.1
        self.n = int(round(self.l/self.dl))
        dummyList = self.n*[1.0]
        self.curvature = array.array("d")
        self.curvature.extend(dummyList)

    # TODO: should this be updated from within the class? e.g. call the update every time-step? Probably not

@dataclass
class BodyLinearVelocity(PythonMsg):
    v_long: float = field(default=0)
    v_tran: float = field(default=0)
    v_n: float = field(default=0)

@dataclass
class BodyAngularVelocity(PythonMsg):
    w_phi: float = field(default=0)
    w_theta: float = field(default=0)
    w_psi: float = field(default=0)

@dataclass
class BodyLinearAcceleration(PythonMsg):
    a_long: float = field(default=0)
    a_tran: float = field(default=0)
    a_n: float = field(default=0)

@dataclass
class BodyAngularAcceleration(PythonMsg):
    a_phi: float = field(default=0)
    a_theta: float = field(default=0)
    a_psi: float = field(default=0)

@dataclass
class OrientationEuler(PythonMsg):
    phi: float = field(default=0)
    theta: float = field(default=0)
    psi: float = field(default=0)

@dataclass
class OrientationQuaternion(PythonMsg):
    qr: float = field(default=1)
    qi: float = field(default=0)
    qj: float = field(default=0)
    qk: float = field(default=0)

    def e1(self):
        return np.array([[1 - 2 * self.qj ** 2 - 2 * self.qk ** 2,
                          2 * (self.qi * self.qj + self.qk * self.qr),
                          2 * (self.qi * self.qk - self.qj * self.qr)]]).T

    def e2(self):
        return np.array([[2 * (self.qi * self.qj - self.qk * self.qr),
                          1 - 2 * self.qi ** 2 - 2 * self.qk ** 2,
                          2 * (self.qj * self.qk + self.qi * self.qr)]]).T

    def e3(self):
        return np.array([[2 * (self.qi * self.qk + self.qj * self.qr),
                          2 * (self.qj * self.qk - self.qi * self.qr),
                          1 - 2 * self.qi ** 2 - 2 * self.qj ** 2]]).T

    def R(self):
        return np.array([[1 - 2 * self.qj ** 2 - 2 * self.qk ** 2, 2 * (self.qi * self.qj - self.qk * self.qr),
                          2 * (self.qi * self.qk + self.qj * self.qr)],
                         [2 * (self.qi * self.qj + self.qk * self.qr), 1 - 2 * self.qi ** 2 - 2 * self.qk ** 2,
                          2 * (self.qj * self.qk - self.qi * self.qr)],
                         [2 * (self.qi * self.qk - self.qj * self.qr), 2 * (self.qj * self.qk + self.qi * self.qr),
                          1 - 2 * self.qi ** 2 - 2 * self.qj ** 2]])

    def Rinv(self):
        return np.array([[1 - 2 * self.qj ** 2 - 2 * self.qk ** 2, 2 * (self.qi * self.qj + self.qk * self.qr),
                          2 * (self.qi * self.qk - self.qj * self.qr)],
                         [2 * (self.qi * self.qj - self.qk * self.qr), 1 - 2 * self.qi ** 2 - 2 * self.qk ** 2,
                          2 * (self.qj * self.qk + self.qi * self.qr)],
                         [2 * (self.qi * self.qk + self.qj * self.qr), 2 * (self.qj * self.qk - self.qi * self.qr),
                          1 - 2 * self.qi ** 2 - 2 * self.qj ** 2]])

    def norm(self):
        return np.sqrt(self.qr ** 2 + self.qi ** 2 + self.qj ** 2 + self.qk ** 2)

    def normalize(self):
        norm = self.norm()

        self.qr /= norm
        self.qi /= norm
        self.qj /= norm
        self.qk /= norm
        return

    def from_yaw(self, yaw):
        self.qi = 0
        self.qj = 0
        self.qr = np.cos(yaw / 2)
        self.qk = np.sin(yaw / 2)
        self.normalize()

        if abs(yaw - self.to_yaw()) % (2 * np.pi) > 1e-9: pdb.set_trace()

        return

    def to_yaw(self):
        return 2 * np.arctan2(self.qk, self.qr)

    def qdot(self, w: BodyAngularVelocity) -> 'OrientationQuaternion':
        qdot = OrientationQuaternion()
        qdot.qr = -0.5 * (self.qi * w.w_phi + self.qj * w.w_theta + self.qk * w.w_psi)
        qdot.qi = 0.5 * (self.qr * w.w_phi + self.qj * w.w_psi - self.qk * w.w_theta)
        qdot.qj = 0.5 * (self.qr * w.w_theta + self.qk * w.w_phi - self.qi * w.w_psi)
        qdot.qk = 0.5 * (self.qr * w.w_psi + self.qi * w.w_theta - self.qj * w.w_phi)
        return qdot

@dataclass
class ParametricPose(PythonMsg):
    s: float = field(default = 0)
    x_tran: float = field(default = 0)
    n: float = field(default = 0)
    e_psi: float = field(default = 0)

@dataclass
class ParametricVelocity(PythonMsg):
    ds: float = field(default = 0)
    dx_tran: float = field(default = 0)
    dn: float = field(default = 0)
    de_psi: float = field(default = 0)

@dataclass
class VehicleState(PythonMsg):
    '''
    Complete vehicle state (local, global, and input)
    '''
    t: float  = field(default = None)    # time in seconds

    x: Position = field(default=None)  # global position

    v: BodyLinearVelocity = field(default=None)  # body linear velocity
    w: BodyAngularVelocity = field(default=None)  # body angular velocity
    a: BodyLinearAcceleration = field(default=None)  # body linear acceleration
    aa: BodyAngularAcceleration = field(default=None)  # body angular acceleration

    q: OrientationQuaternion = field(default=None)  # global orientation (qr, qi, qj, qk)
    e: OrientationEuler = field(default=None)       # global orientation (phi, theta, psi)

    p: ParametricPose = field(default=None)  # parametric position (s,y, ths)
    pt: ParametricVelocity = field(default=None)  # parametric velocity (ds, dy, dths)

    u: VehicleActuation = field(default=None)

    lookahead: TrackLookahead = field(default=None) # TODO Find a good field name :(

    lap_num: int = field(default = None)

    #TODO: Maybe need a different way of storing the covariance matricies
    # For covariances, only store the upper triangular part
    local_state_covariance: array.array = field(default = None) # Vectorized upper triangular part of covariance matrix with main diagonal order [v_long, v_tran, psidot, e_psi, s, e_y]
    global_state_covariance: array.array = field(default = None) # Vectorized upper triangular part of covariance matrix with main diagonal order [x, y, psi, v_long, v_tran, psidot]

    def __post_init__(self):
        if self.x is None: self.x = Position()
        if self.u is None: self.u = VehicleActuation()
        if self.lookahead is None: self.lookahead = TrackLookahead()
        if self.v is None: self.v = BodyLinearVelocity()
        if self.w is None: self.w = BodyAngularVelocity()
        if self.a is None: self.a = BodyLinearAcceleration()
        if self.aa is None: self.aa = BodyAngularAcceleration()
        if self.q is None: self.q = OrientationQuaternion()
        if self.e is None: self.e = OrientationEuler()
        if self.p is None: self.p = ParametricPose()
        if self.pt is None: self.pt = ParametricVelocity()
        return

    def get_R(self, reverse = False):
        #Warning - Not suitable for general 3D case
        psi = self.psi
        return np.array([[np.cos(psi), -np.sin(psi),0],[np.sin(psi),np.cos(psi),0],[0,0,0]])

    # Single operation with rotation matrix
    # def update_body_velocity_from_global(self):
    #     self.v_long =  self.v_x * np.cos(self.psi) + self.v_y * np.sin(self.psi)
    #     self.v_tran = -self.v_x * np.sin(self.psi) + self.v_y * np.cos(self.psi)
    #     self.a_long =  self.a_x * np.cos(self.psi) + self.a_y * np.sin(self.psi)
    #     self.a_tran = -self.a_y * np.sin(self.psi) + self.a_y + np.cos(self.psi)
    #
    # def update_global_velocity_from_body(self):
    #     self.v_x =  self.v_long * np.cos(self.psi) - self.v_tran * np.sin(self.psi)
    #     self.v_y =  self.v_long * np.sin(self.psi) + self.v_tran * np.cos(self.psi)
    #     self.a_x =  self.a_long * np.cos(self.psi) - self.a_tran * np.sin(self.psi)
    #     self.a_y =  self.a_long * np.sin(self.psi) + self.a_tran * np.cos(self.psi)

    def pack_list(self, use_numpy = False):
        '''
        Takes in a list or array of VehicleState objects and creates a single VehicleState object where each field is a np.array of the fields of the original list

        usage: pytypes.VehicleState.pack_list(array_of_VehicleState)
        '''
        msg = VehicleState()
        N = len(self)
        for key in vars(msg).keys():
            msg.__setattr__(key,[])

            for i in range(N):
                msg.__getattribute__(key).append(self[i].__getattribute__(key))
            if use_numpy:
                msg.__setattr__(key, np.array(msg.__getattribute__(key)))

        return msg

    def copy_control(self, destination):
        '''
        copies control state form self to destination
        '''
        destination.t = self.t
        destination.u_a = self.u.u_a
        destination.u_steer = self.u.u_steer
        return

# TODO: Change to array of VehicleState
@dataclass
class VehiclePrediction(PythonMsg):
    '''
    Complete vehicle coordinates (local, global, and input)
    '''
    t: float  = field(default = None)    # time in seconds

    x: array.array = field(default = None)     # global x coordinate in meters
    y: array.array = field(default = None)     # global y coordinate in meters

    v_x: array.array = field(default = None)   # global x velocity in m/s
    v_y: array.array = field(default = None)   # global y velocity in m/s

    a_x: array.array = field(default = None)   # global x acceleration in m/s^2
    a_y: array.array = field(default = None)   # global y acceleration in m/s^2

    psi: array.array = field(default = None)      # global vehicle heading angle
    psidot: array.array = field(default = None)   # global and local angular velocity of car

    v_long: array.array = field(default = None)   # longitudinal velocity (in the direction of psi)
    v_tran: array.array = field(default = None)   # transverse velocity   (orthogonal to the direction of psi)

    a_long: array.array = field(default = None)   # longitudinal velocity (in the direction of psi)
    a_tran: array.array = field(default = None)   # transverse velocity   (orthogonal to the direction of psi)

    e_psi: array.array = field(default = None)    # heading error between car and track
    s: array.array = field(default = None)        # path length along center of track to projected position of car
    x_tran: array.array = field(default = None)   # deviation from centerline (transverse position)

    u_a: array.array = field(default = None)      # acceleration output
    u_steer: array.array = field(default = None)  # steering angle output

    lap_num: int = field(default = None)

    #TODO: Maybe need a different way of storing the covariance matricies
    # For covariances, only store the upper triangular part
    local_state_covariance: array.array = field(default = None) # Vectorized upper triangular part of covariance matrix with main diagonal order [v_long, v_tran, psidot, e_psi, s, e_y]
    global_state_covariance: array.array = field(default = None) # Vectorized upper triangular part of covariance matrix with main diagonal order [x, y, psi, v_long, v_tran, psidot]

    def update_body_velocity_from_global(self):
        self.v_long =  (np.multiply(self.v_x, np.cos(self.psi)) + np.multiply(self.v_y, np.sin(self.psi))).tolist()
        self.v_tran = (-np.multiply(self.v_x, np.sin(self.psi)) + np.multiply(self.v_y, np.cos(self.psi))).tolist()
        self.a_long =  (np.multiply(self.a_x, np.cos(self.psi)) + np.multiply(self.a_y, np.sin(self.psi))).tolist()
        self.a_tran = (-np.multiply(self.a_y, np.sin(self.psi)) + np.multiply(self.a_y, np.cos(self.psi))).tolist()

    def update_global_velocity_from_body(self):
        self.v_x =  (np.multiply(self.v_long, np.cos(self.psi)) - np.multiply(self.v_tran, np.sin(self.psi))).tolist()
        self.v_y =  (np.multiply(self.v_long, np.sin(self.psi)) + np.multiply(self.v_tran, np.cos(self.psi))).tolist()
        self.a_x =  (np.multiply(self.a_long, np.cos(self.psi)) - np.multiply(self.a_tran, np.sin(self.psi))).tolist()
        self.a_y =  (np.multiply(self.a_long, np.sin(self.psi)) + np.multiply(self.a_tran, np.cos(self.psi))).tolist()

max_stars_repo_name: zhu-edward/DGSQP

max_issues_repo_name: zhu-edward/DGSQP

max_forks_repo_name: zhu-edward/DGSQP


text:
import numpy as np
import glob,os
import pandas as pd
import matplotlib as mpl
mpl.use('Agg')
from matplotlib import pyplot as plt
from matplotlib import cm
import seaborn as sns

a = 10**np.linspace(-1,2,10)
b = 10**np.linspace(-2,1,10)
x = np.zeros((10,10))
FP = glob.glob(os.path.join('Logs/Optimize.o19226927.*'))
for fp in FP:
	f = open(fp)
	L = f.readlines()
	f.close()
	if len(L) > 0:
		i = [float(l.strip().split('\t')[-1]) for l in L if 'alpha_cond' in l][0]
		j = [float(l.strip().split('\t')[-1]) for l in L if 'alpha_corr' in l][0]
		x[np.where(abs(a-i) < 1e-4)[0][0],np.where(abs(b-j) < 1e-4)[0][0]] = float(L[-2].strip().split()[-1])

astr = ['%.2f'%i for i in a]
df = pd.DataFrame(x,columns=b,index=astr)
sns_plt = sns.heatmap(df,mask=(df == 0),annot=True)
plt.savefig('alpha_grid.png')
plt.close()

max_stars_repo_name: cleary-lab/CISI

max_issues_repo_name: cleary-lab/CISI

max_forks_repo_name: cleary-lab/CISI