text:
import os
import pathlib as path
import subprocess
import tempfile

import numpy as np

from . import helpers

class Species:
    """
    A single bi-Maxwellian species.

    Parameters
    ----------
    q : float
        Particle charge, in units of proton charge. e.g. for a proton
        enter ``1``, for an electron enter ``-1``.
    m : float
        Particle mass, in units of proton mass. e.g. for a proton
        enter ``1``.
    n : float
        Number density, relative to the number density used to define the
        Alfvén speed.
    v_d : float
        Drift speed as a fraction of the Alfvén speed.
    t_ani : float
        Temperature anisotropy (perpendicular temperature divided by
        parallel temperature).
    beta_par : float
        Parallel beta (thermal pressure divided by magnetic pressure).
    """
    def __init__(self, q, m, n, v_d, t_ani, beta_par):
        assert m > 0, 'Mass must be greater than zero'
        self.q = q
        self.m = m
        self.n = n
        self.v_d = v_d
        self.t_ani = t_ani
        self.beta_par = beta_par

    def __repr__(self):
        return ('Species object with ' +
                'q={self.q}, m={self.m}, n={self.n}, '.format(self=self) +
                'v_d={self.v_d}, t_ani={self.t_ani}, beta_par={self.beta_par}'.format(self=self)
                )

class InputParams:
    """
    Input parameters for calculating a dispersion relation.

    Parameters
    ----------
    species : list of :class:`Species`
        List of the particle species for which to calculate the dispersion.
        **Note** that the order of the species matters, as some of the
        parameters below are defined relative to the first species in the
        *species* list.
    omega_guess : complex
        Initial guess for the frequency, normalised to the gyro frequency of
        the first species.
    propagation_angle : float
        Propagation angle to calculate dispersion relation for. Defined
        with respect to the magnetic field.
    va : float
        Ratio of the Alfvén speed to the speed of light. The number density
        in the Alfvén speed is taken from the first species.
    kzmin : float
        Start of kz range.
    kzmax : float
        End of kz range.

    Other parameters
    ----------------
    kzsteps : int, optional
        Number of steps in kz to calculate dispersion relation at. Points are
        linearly spaced between *kzmin* and *kzmax*.
    numiter : int, optional
        Maximum number of iterations in the Newton method
    det_D_threshold : float, optional
        Threshold for the determinant of the dispersion tensor.
        If det D <= det_D_threshold, the Newton iteration will be stopped.
    n_bessel : int, optional
        Maximum of sum in Bessel functions (both regular and modified).
        Can be very low (e.g., 3) for quasi-parallel propagation.
    bessel_zero : float, optional
        If I_n is less than this value, higher n are neglected.
    vxsteps, vysteps, vzsteps : int, optional
        Steps in the x,y,z directions.
    """
    def __init__(self, species, omega_guess, propagation_angle, va,
                 kzmin, kzmax, kzsteps=200,
                 numiter=1000, det_D_threshold=1e-16,
                 n_bessel=1000, bessel_zero=1e-50,
                 vxsteps=100, vysteps=100, vzsteps=100
                 ):
        self.species = species
        self.omega_guess = omega_guess
        self.propagation_angle = propagation_angle
        self.va = va
        self.kzmin = kzmin
        self.kzmax = kzmax
        self.kzsteps = kzsteps
        self.numiter = numiter
        self.det_D_threshold = det_D_threshold
        self.n_bessel = n_bessel
        self.bessel_zero = bessel_zero
        self.vxsteps = vxsteps
        self.vysteps = vysteps
        self.vzsteps = vzsteps

    @property
    def nspecies(self):
        '''
        Number of species.
        '''
        return len(self.species)

    @property
    def total_charge(self):
        """
        Total charge density. Result is normalised to the first species' charge
        and number density.
        """
        q = 0
        for species in self.species:
            q += species.q * species.n
        return q

    @property
    def total_current(self):
        """
        Total current density.Result is normalised to the first species'
        charge, number density, and the Alfvén speed.
        """
        j = 0
        for species in self.species:
            j += species.q * species.n * species.v_d
        return j

class Result:
    '''
    Result of running the dispersion solver.
    '''
    def __init__(self, input, run_output):
        self.input = input
        self.output = run_output

    @property
    def kz(self):
        '''
        Wave vector normalised to the proton intertial length (:math:`kd_{p}`).
        '''
        return self.output[:, 0]

    @property
    def omega_real(self):
        r'''
        Real part of frequency normalised to the proton gyro-frequency
        (:math:`\omega / \Omega_{p}`).
        '''
        return self.output[:, 2]

    @property
    def omega_imag(self):
        r'''
        Imaginary part of frequency normalised to the proton gyro-frequency
        (:math:`\gamma / \Omega_{p}`).
        '''
        return self.output[:, 3]

    @property
    def EyEx(self):
        '''
        The electric field component ratio :math:`E_{y} / E_{x}`.
        Returned as complex numbers.
        '''
        return self.output[:, 4] + self.output[:, 5] * 1j

    @property
    def EzEx(self):
        '''
        The electric field component ratio :math:`E_{z} / E_{x}`.
        Returned as complex numbers.
        '''
        return self.output[:, 6] + self.output[:, 7] * 1j

def format_input_file(input):
    '''
    Function to create input file.

    Parameters
    ----------
    input : InputParams
    '''
    def create_species_list(attr):
        '''
        Function to populate species properties.
        '''
        n = 0
        ret = ''
        for species in input.species:
            x = getattr(species, attr)
            ret += '{},'.format(helpers._float_to_fortran_str(x))
            n += 1
        for i in range(n, 10):
            ret += '0.d0,'
        ret = ret[:-1]
        return ret

    input_file_sting = r'''
&parameters

! Number of species
! (up to 10 species are possible)
numspec={numspec}

! Maximum number of iterations in the Newton method
numiter={numiter}

! Threshold for the determinant of the dispersion tensor:
! If det D <= det_D_threshold, the Newton iteration will be stopped
det_D_threshold={detD}

! Maximum of sum in Bessel functions (both regular and modified)
! can be very low (e.g., 3) for quasi-parallel propagation
nmax={nmax}

! If I_n is less than this value, higher n are neglected:
Bessel_zero={bessel_zero}

! Initial guess for the frequency (complex number):
initial_guess=({omega_guess_real},{omega_guess_imag})

! Range of values to scan over in kz:
krange={kzmin},{kzmax}

! Number of steps to scan over kzrange:
ksteps={kzsteps}

! Temperature anisotropy (Tperp/Tparallel)
alpha={alpha}

! Parallel beta of the species
beta={beta}

! Charge of the species in units of the first ion charge
charge={charge}

! Mass of the species in units of ion mass
mass={mass}

! Density of the species in units of proton density
density={density}

! Drift speed of the species in units of proton Alfven speed
vdrift={vdrift}

! Range of angle of propagation (in degrees)
theta_range={propagation_angle}, {propagation_angle}

! Steps in the theta direction
theta_steps=1

! Alfven speed divided by speed of light
vAc={va}

! Amplitude mode: If 1, then ampl = delta B_x / B0. If 2, then ampl = delta B_y / B0.
! If 3, then ampl = delta B_z / B0
ampl_mode=3

! Amplitude for the calculation of polarization properties:
ampl=1.d0

! Print warnings to std. output and stop program in case of problems
output_warning=.FALSE.

!! The following parameters are used to determine delta f in write_delta_f:
! Maximum of sum in Bessel function:
! can be very low (e.g., 3) for quasi-parallel propagation
mmax=1000

! If I_n is less than this value, higher m are neglected:
Bessel_zero_deltaf=1.d-50

! Steps in vx, vy, and vz (standard: 100):
vxsteps={vxsteps}
vysteps={vysteps}
vzsteps={vzsteps}

! Range in vpar and vperp:
vxrange=-1.d0,1.d0
vyrange=-1.d0,1.d0
vzrange=-1.d0,1.d0

! Number of time steps for one full period (standard: 25):
timesteps=40

! If periods is .TRUE., then num_periods is number of periods (2 PI / omega_r).
! If periods is .FALSE., then num_periods is number of gyro-periods (2 PI / Omega_p).
periods=T

! Number of periods (standard: 8):
num_periods=1

! Include wave damping in delta f calculation or not (exp(gamma*t)-term):
damping=F

! If const_r is .TRUE., then delta f is evaluated at point r = 0 -> cos(-omega * t)
! If const_r is .FALSE., then the wave period is a function of space at -> cos(k * r)
!		with k * r from 0 ... 2 pi at a fixed time. We "fly" along the k-direction over one wave train.
! 		If set on .FALSE., periods does not make a difference
const_r=T

! Write moments to a file or not:
output_mom=F

! Write E & B vectors to a file or not:
output_EB=F

! Override the k,theta values specified above with the ones read from a file:
kth_file=F

! Filename to read k, theta values from:
! kth_filename='kth_kawtc'
/
    '''.format(numspec=input.nspecies,
               numiter=input.numiter,
               detD=helpers._float_to_fortran_str(input.det_D_threshold),
               nmax=input.n_bessel,
               bessel_zero=helpers._float_to_fortran_str(input.bessel_zero),
               omega_guess_real=helpers._float_to_fortran_str(input.omega_guess.real),
               omega_guess_imag=helpers._float_to_fortran_str(input.omega_guess.imag),
               propagation_angle=helpers._float_to_fortran_str(input.propagation_angle),
               va=helpers._float_to_fortran_str(input.va),
               vxsteps=input.vxsteps,
               vysteps=input.vysteps,
               vzsteps=input.vzsteps,
               alpha=create_species_list('t_ani'),
               beta=create_species_list('beta_par'),
               charge=create_species_list('q'),
               mass=create_species_list('m'),
               density=create_species_list('n'),
               vdrift=create_species_list('v_d'),
               kzmin=helpers._float_to_fortran_str(input.kzmin),
               kzmax=helpers._float_to_fortran_str(input.kzmax),
               kzsteps=input.kzsteps
               )
    return input_file_sting

def run(input):
    '''
    Run the dispersion solver for a given input.

    Parameters
    ----------
    input : InputParams
    '''
    abspath = path.Path(os.path.dirname(__file__))
    nhds_folder = abspath / 'NHDS'
    (nhds_folder / 'obj').mkdir(exist_ok=True)
    (nhds_folder / 'bin').mkdir(exist_ok=True)
    parameter_file = nhds_folder / 'parameters.in'
    # Create input file
    with open(parameter_file, mode='w') as f:
        f.write(format_input_file(input))
    # Check if NHDS is compiled
    binary = nhds_folder / 'NHDS'
    if not binary.exists():
        helpers.compile_nhds(nhds_folder)
    subprocess.check_output(
        ['./NHDS', 'parameters.in'], universal_newlines=True, cwd=nhds_folder)

    output_file = nhds_folder / 'output_parameters.in.dat'
    out = np.loadtxt(output_file, skiprows=1)
    return Result(input, out)

from ._version import get_versions
__version__ = get_versions()['version']
del get_versions

max_stars_repo_name: dstansby/nhdspy

max_issues_repo_name: dstansby/NHDSPy

max_forks_repo_name: dstansby/nhdspy


text:
#!/usr/bin/env python
# Copyright (c) 2016 The UUV Simulator Authors.
# All rights reserved.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

import math
import numpy
import rospy
import tf.transformations as trans
from PID import PIDRegulator

from dynamic_reconfigure.server import Server
from uuv_control_cascaded_pid.cfg import PositionControlConfig
import geometry_msgs.msg as geometry_msgs
from nav_msgs.msg import Odometry
from rospy.numpy_msg import numpy_msg

class PositionControllerNode:
    def __init__(self):
        print('PositionControllerNode: initializing node')

        self.config = {}

        self.pos_des = numpy.zeros(3)
        self.quat_des = numpy.array([0, 0, 0, 1])

        self.initialized = False

        # Initialize pids with default parameters
        self.pid_rot = PIDRegulator(1, 0, 0, 1)
        self.pid_pos = PIDRegulator(1, 0, 0, 1)

        # ROS infrastructure
        self.sub_cmd_pose = rospy.Subscriber('cmd_pose', numpy_msg(geometry_msgs.Pose), self.cmd_pose_callback)
        self.sub_odometry = rospy.Subscriber('odom', numpy_msg(Odometry), self.odometry_callback)
        self.pub_cmd_vel = rospy.Publisher('cmd_vel', geometry_msgs.Twist, queue_size=10)
        self.srv_reconfigure = Server(PositionControlConfig, self.config_callback)

    def cmd_pose_callback(self, msg):
        """Handle updated set pose callback."""
        # Just store the desired pose. The actual control runs on odometry callbacks
        p = msg.position
        q = msg.orientation
        self.pos_des = numpy.array([p.x, p.y, p.z])
        self.quat_des = numpy.array([q.x, q.y, q.z, q.w])

    def odometry_callback(self, msg):
        """Handle updated measured velocity callback."""

        if not bool(self.config):
            return

        p = msg.pose.pose.position
        q = msg.pose.pose.orientation
        p = numpy.array([p.x, p.y, p.z])
        q = numpy.array([q.x, q.y, q.z, q.w])

        if not self.initialized:
            # If this is the first callback: Store and hold latest pose.
            self.pos_des  = p
            self.quat_des = q
            self.initialized = True

        # Compute control output:
        t = msg.header.stamp.to_sec()

        # Position error
        e_pos_world = self.pos_des - p
        e_pos_body = trans.quaternion_matrix(q).transpose()[0:3,0:3].dot(e_pos_world)

        # Error quaternion wrt body frame
        e_rot_quat = trans.quaternion_multiply(trans.quaternion_conjugate(q), self.quat_des)

        if numpy.linalg.norm(e_pos_world[0:2]) > 5.0:
            # special case if we are far away from goal:
            # ignore desired heading, look towards goal position
            heading = math.atan2(e_pos_world[1],e_pos_world[0])
            quat_des = numpy.array([0, 0, math.sin(0.5*heading), math.cos(0.5*heading)])
            e_rot_quat = trans.quaternion_multiply(trans.quaternion_conjugate(q), quat_des)

        # Error angles
        e_rot = numpy.array(trans.euler_from_quaternion(e_rot_quat))

        v_linear = self.pid_pos.regulate(e_pos_body, t)
        v_angular = self.pid_rot.regulate(e_rot, t)

        # Convert and publish vel. command:
        cmd_vel = geometry_msgs.Twist()
        cmd_vel.linear = geometry_msgs.Vector3(*v_linear)
        cmd_vel.angular = geometry_msgs.Vector3(*v_angular)
        self.pub_cmd_vel.publish(cmd_vel)

    def config_callback(self, config, level):
        """Handle updated configuration values."""
        # Config has changed, reset PID controllers
        self.pid_pos = PIDRegulator(config['pos_p'], config['pos_i'], config['pos_d'], config['pos_sat'])
        self.pid_rot = PIDRegulator(config['rot_p'], config['rot_i'], config['rot_d'], config['rot_sat'])

        self.config = config

        return config

if __name__ == '__main__':
    print('starting PositionControl.py')
    rospy.init_node('position_control')

    try:
        node = PositionControllerNode()
        rospy.spin()
    except rospy.ROSInterruptException:
        print('caught exception')
    print('exiting')

ext: py

max_stars_repo_name: shubhamkorde/AnahitaPlus

max_issues_repo_name: shubhamkorde/AnahitaPlus

max_forks_repo_name: shubhamkorde/AnahitaPlus

reason: import numpy


text:
import re
import pandas as pd
import numpy as np
import sklearn
from sklearn import preprocessing
from sklearn.preprocessing import OneHotEncoder, LabelEncoder
from sklearn.model_selection import train_test_split
from sklearn.feature_extraction.text import TfidfVectorizer, CountVectorizer
import nltk
from nltk.corpus import stopwords
from nltk.tokenize import word_tokenize
from nltk.stem import WordNetLemmatizer
import mlflow
import mlflow.sklearn
from sklearn.pipeline import Pipeline

df =  pd.read_csv('ISEAR.csv', names=['emotion','content'])

#define le as LabelEncoder to encode the categories
le = preprocessing.LabelEncoder()

#encode emotion
df['emotion'] = le.fit_transform(df['emotion'])

#print labels assigened to emotion
list(le.inverse_transform([0,1, 2, 3]))

#import libraries and initiate required objects for processing
tokenizer = nltk.tokenize.TreebankWordTokenizer()
stop_words = stopwords.words("english")
from nltk.stem.snowball import SnowballStemmer
stemmer = SnowballStemmer("english")
lemmatizer = WordNetLemmatizer()
def process_text(text):
    #convert text to lowercase
    text= text.apply(lambda x:x.lower())
    #remove multiple spaces from text
    text= text.apply(lambda x:re.sub(' +', ' ',x))
    # remove alpha numeric characeter from text using regex
    text= text.apply(lambda x:x.lower())
    text.str.replace(r"\d+", "")
    text.str.replace('[^\w\s]','')
    text.str.replace(r"[︰-＠]", "")
    text.str.replace(r"", "")
    text.str.replace('\d+', '')
    #tokenize
    text = text.apply(tokenizer.tokenize)
    #remove english stop words form text
    text =text.apply(lambda x: [item for item in x if item not in stop_words])
    #stemming text
    text =text.apply(lambda x: [stemmer.stem(e) for e in x])
    #stripping spaces from items of array of texts
    text=text.apply(lambda x: [e.strip() for e in x])
    #remove all the characters from array of texts if the length of item is 1
    text=text.apply( lambda x: [ y for y in x if len(y)>1 ])
    #remove all the digits from the text
    text=text.apply( lambda x: [ y for y in x if not y.isdigit()] )
    # lemattizing text
    text = text.apply(lambda x: ' '.join(lemmatizer.lemmatize(token) for token in x))
    return text

    df['content'] = process_text(df['content'])


text:
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Jul  6 15:55:12 2021

@author: mike_ubuntu
"""
import numpy as np
from typing import Union, Callable
from collections import OrderedDict
import warnings

import epde.globals as global_var

from epde.moeadd.moeadd import *
from epde.moeadd.moeadd_supplementary import *

from epde.prep.DomainPruning import Domain_Pruner

from epde.operators.ea_stop_criteria import Stop_condition, Iteration_limit
import epde.operators.sys_search_operators as operators

from epde.evaluators import simple_function_evaluator, trigonometric_evaluator
from epde.supplementary import Define_Derivatives
from epde.cache.cache import upload_simple_tokens, upload_grids, prepare_var_tensor, np_ndarray_section
from epde.prep.derivatives import Preprocess_derivatives
from epde.eq_search_strategy import Strategy_director, Strategy_director_solver
from epde.structure import Equation

from epde.interface.token_family import TF_Pool, Token_family
from epde.interface.type_checks import *
from epde.interface.prepared_tokens import Prepared_tokens, Custom_tokens

class Input_data_entry(object):
    def __init__(self, var_name, data_tensor, coord_tensors = None):
        self.var_name = var_name
        check_nparray(data_tensor)
        if coord_tensors is not None:
            check_nparray_iterable(coord_tensors)
            if any([tensor.shape != data_tensor.shape for tensor in coord_tensors]):
                print('data_tensor.shape', data_tensor.shape, 'tensors:',
                      [tensor.shape for tensor in coord_tensors])
                raise ValueError('mismatching shapes of coordinate tensors and input data')
            self.coord_tensors = coord_tensors
        else:
            axes = []
            for ax_idx in range(data_tensor.ndim):
                axes.append(np.linspace(0., 1., data_tensor.shape[ax_idx]))
            self.coord_tensors = np.meshgrid(*axes)
        self.data_tensor = data_tensor

    # def set_derivatives(self, deriv_tensors = None, data_filename = None,
    #                     deriv_filename = None, smooth = True, sigma = 5, max_order = 1, mp_poolsize = 2):
#        coord_names = ['x' + str(coord_idx) for coord_idx in range(len(self.coord_tensors))]
    def set_derivatives(self, deriv_tensors = None, method = 'ANN', max_order = 1, method_kwargs = {}):

        deriv_names, deriv_orders = Define_Derivatives(self.var_name, dimensionality=self.data_tensor.ndim,
                                                       max_order = max_order)

        self.names = deriv_names # coord_names +  Define_Derivatives(self.var_name, dimensionality=self.data_tensor.ndim,
                                    #                  max_order = max_order)
        self.d_orders = deriv_orders
        if deriv_tensors is None:
            method_kwargs['max_order'] = max_order
            if self.coord_tensors is not None and 'grid' not in method_kwargs.keys():
                method_kwargs['grid'] = self.coord_tensors
            _, self.derivatives = Preprocess_derivatives(self.data_tensor, method=method,
                                                         method_kwargs=method_kwargs)
            self.deriv_properties = {'max order' : max_order,
                                     'dimensionality' : self.data_tensor.ndim}
        else:
            self.derivatives = deriv_tensors
            self.deriv_properties = {'max order' : max_order,
                                     'dimensionality' : self.data_tensor.ndim}

    def use_global_cache(self, grids_as_tokens = True, set_grids = True, memory_for_cache = 5,
                         boundary : int = 0):

        print(type(self.data_tensor), type(self.derivatives))
        derivs_stacked = prepare_var_tensor(self.data_tensor, self.derivatives,
                                            time_axis = global_var.time_axis, boundary = boundary)
#                                            axes = self.coord_tensors)
        if isinstance(self.coord_tensors, (list, tuple)):
            coord_tensors_cut = []
            for tensor in self.coord_tensors:
#                print('appending grid of shape', tensor.shape, 'cut to ', np_ndarray_section(tensor, boundary = boundary))
                coord_tensors_cut.append(np_ndarray_section(tensor, boundary = boundary))
        elif isinstance(self.coord_tensors, np.ndarray):
            coord_tensors_cut = np_ndarray_section(self.coord_tensors, boundary = boundary)
        else:
            raise TypeError('Coordinate tensors are presented in format, other than np.ndarray or list/tuple of np.ndarray`s')

        try:
#            print(self.names, derivs_stacked.shape)
            upload_simple_tokens(self.names, global_var.tensor_cache, derivs_stacked)
            if set_grids:
    #                print('setting grids')
                memory_for_cache = int(memory_for_cache/2)
    #                global_var.grid_cache.memory_usage_properties(obj_test_case = self.data_tensor,
    #                                                                mem_for_cache_frac = memory_for_cache)
                upload_grids(coord_tensors_cut, global_var.grid_cache)
                print(f'completed grid cache with {len(global_var.grid_cache.memory_default)} tensors with labels {global_var.grid_cache.memory_default.keys()}')

        except AttributeError:
            print('Somehow, cache had not been initialized')
            print(self.names, derivs_stacked.shape)
            print(global_var.tensor_cache.memory_default.keys())
            global_var.init_caches(set_grids = set_grids)
            if set_grids:
                print('setting grids')
                memory_for_cache = int(memory_for_cache/2)
                global_var.grid_cache.memory_usage_properties(obj_test_case = self.data_tensor,
                                                                mem_for_cache_frac = memory_for_cache)
                upload_grids(coord_tensors_cut, global_var.grid_cache)
                print(f'completed grid cache with {len(global_var.grid_cache.memory_default)} tensors with labels {global_var.grid_cache.memory_default.keys()}')
#                print()

            global_var.tensor_cache.memory_usage_properties(obj_test_case = self.data_tensor,
                                                            mem_for_cache_frac = memory_for_cache)
            print(self.names, derivs_stacked.shape)
            upload_simple_tokens(self.names, global_var.tensor_cache, derivs_stacked)

        global_var.tensor_cache.use_structural()

def simple_selector(sorted_neighbors, number_of_neighbors = 4):
    return sorted_neighbors[:number_of_neighbors]

class epde_search(object):
    def __init__(self, use_default_strategy : bool = True, eq_search_stop_criterion : Stop_condition = Iteration_limit,
                 director = None, equation_type : set = {'PDE', 'derivatives only'}, time_axis : int = 0,
                 init_cache : bool = True, example_tensor_shape : Union[tuple, list] = (1000,),
                 set_grids : bool = True, eq_search_iter : int = 300, use_solver : bool = False,
                 dimensionality : int = 1):
        '''

        Intialization of the epde search object. Here, the user can declare the properties of the
        search mechainsm by defining evolutionary search strategy.

        Parameters:
        --------------

        use_default_strategy : bool, optional
            True (base and recommended value), if the default evolutionary strategy shall be used,
            False if the user - defined strategy will be passed further.

        eq_search_stop_criterion : epde.operators.ea_stop_criteria.Stop_condition object, optional
            The stop condition for the evolutionary search of the equation. Default value
            represents loop over 300 evolutionary epochs.

        director : obj, optional
            User-defined director, responsible for construction of evolutionary strategy;
            shall be declared for very specific tasks.

        equation_type : set of str, optional
            Will define equation type, TBD later.

        '''
        global_var.set_time_axis(time_axis)
#        eq_search_stop_criterion_obj = eq_search_stop_criterion()
        if director is not None and not use_default_strategy:
            self.director = director
        elif director is None and use_default_strategy:
            if use_solver:
                self.director = Strategy_director_solver(eq_search_stop_criterion, {'limit' : eq_search_iter},
                                                         dimensionality=dimensionality)
            else:
                self.director = Strategy_director(eq_search_stop_criterion, {'limit' : eq_search_iter})
            self.director.strategy_assembly()
        else:
            raise NotImplementedError('Wrong arguments passed during the epde search initialization')
        self.set_moeadd_params()
        self.search_conducted = False

        if init_cache:
            global_var.init_caches(set_grids)

    def set_memory_properties(self, example_tensor, mem_for_cache_frac = None, mem_for_cache_abs = None):
        if global_var.grid_cache is not None:
            if mem_for_cache_frac is not None:
                mem_for_cache_frac = int(mem_for_cache_frac/2.)
            else:
                mem_for_cache_abs = int(mem_for_cache_abs/2.)
            global_var.tensor_cache.memory_usage_properties(example_tensor, mem_for_cache_frac, mem_for_cache_abs)
        global_var.tensor_cache.memory_usage_properties(example_tensor, mem_for_cache_frac, mem_for_cache_abs )

    def set_moeadd_params(self, population_size : int = 6, solution_params : dict = {},
                          delta : float = 1/50., neighbors_number : int = 3,
                          NDS_method : Callable = fast_non_dominated_sorting,
                          NDL_update_method : Callable = NDL_update,
                          subregion_mating_limitation : float = .95,
                          PBI_penalty : float = 1., training_epochs : int = 100,
                          neighborhood_selector : Callable = simple_selector,
                          neighborhood_selector_params : tuple = (4,)):
        '''

        Setting the parameters of the multiobjective evolutionary algorithm. declaration of
        the default values is held in the initialization of epde_search object.

        Parameters:
        ------------

        population_size : int, optional
            The size of the population of solutions, created during MO - optimization, default 6.

        solution_params : dict, optional
            Dictionary, containing additional parameters to be sent into the newly created
            solutions.

        delta : float, optional
            parameter of uniform spacing between the weight vectors; *H = 1 / delta*
            should be integer - a number of divisions along an objective coordinate axis.

        neighbors_number : int, *> 0*, optional
            number of neighboring weight vectors to be considered during the operation
            of evolutionary operators as the "neighbors" of the processed sectors.

        NDS_method : function, optional, default ``moeadd.moeadd_supplementary.fast_non_dominated_sorting``
            Method of non-dominated sorting of the candidate solutions. The default method is implemented according to the article
            *K. Deb, A. Pratap, S. Agarwal, and T. Meyarivan, “A fast and elitist
            multiobjective genetic algorithm: NSGA-II,” IEEE Trans. Evol. Comput.,
            vol. 6, no. 2, pp. 182–197, Apr. 2002.*

        NDL_update : function, optional, defalut ``moeadd.moeadd_supplementary.NDL_update``
            Method of adding a new solution point into the objective functions space, introduced
            to minimize the recalculation of the non-dominated levels for the entire population.
            The default method was taken from the *K. Li, K. Deb, Q. Zhang, and S. Kwong, “Efficient non-domination level
            update approach for steady-state evolutionary multiobjective optimization,”
            Dept. Electr. Comput. Eng., Michigan State Univ., East Lansing,
            MI, USA, Tech. Rep. COIN No. 2014014, 2014.*

        neighborhood_selector : function, optional
            Method of finding "close neighbors" of the vector with proximity list.
            The baseline example of the selector, presented in
            ``moeadd.moeadd_stc.simple_selector``, selects n-adjacent ones.

        subregion_mating_limitation : float, optional
            The probability of mating selection to be limited only to the selected
            subregions (adjacent to the weight vector domain). :math:`\delta \in [0., 1.)

        neighborhood_selector_params : tuple/list or None, optional
            Iterable, which will be passed into neighborhood_selector, as
            an arugument. *None*, is no additional arguments are required inside
            the selector.

        epochs : int, optional
            Maximum number of iterations, during that the optimization will be held.
            Note, that if the algorithm converges to a single Pareto frontier,
            the optimization is stopped.

        PBI_penalty :  float, optional
            The penalty parameter, used in penalty based intersection
            calculation, defalut value is 1.

        '''
        self.moeadd_params = {'weights_num' : population_size, 'pop_size' : population_size,
                              'delta' : delta, 'neighbors_number' : neighbors_number,
                              'solution_params': solution_params,
                              'NDS_method' : NDS_method,
                              'NDL_update' : NDL_update_method}

        self.moeadd_optimization_params = {'neighborhood_selector' : neighborhood_selector,
                                           'neighborhood_selector_params' : neighborhood_selector_params,
                                           'delta' : subregion_mating_limitation,
                                           'epochs' : training_epochs,
                                           'PBI_penalty' : PBI_penalty}

    def set_domain_pruning(self, pivotal_tensor_label = None, pruner = None,
                             threshold : float = 1e-5, division_fractions = 3,
                             rectangular : bool = True):
        if pruner is not None:
            self.pruner = pruner
        else:
            self.pruner = Domain_Pruner(domain_selector_kwargs={'threshold' : threshold})
        if not self.pruner.bds_init:
            if pivotal_tensor_label is None:
                pivotal_tensor_label = ('du/dx'+str(global_var.time_axis + 1), (1.0,))
            elif isinstance(pivotal_tensor_label, str):
                pivotal_tensor_label = (pivotal_tensor_label, (1.0,))
            elif not isinstance(pivotal_tensor_label, tuple):
                raise TypeError('Label of the pivotal tensor shall be declared with str or tuple')

            pivotal_tensor = global_var.tensor_cache.get(pivotal_tensor_label)
            self.pruner.get_boundaries(pivotal_tensor, division_fractions = division_fractions,
                                       rectangular = rectangular)

        global_var.tensor_cache.prune_tensors(self.pruner)
        if global_var.grid_cache is not None:
            global_var.grid_cache.prune_tensors(self.pruner)

    def create_pool(self, data : Union[np.ndarray, list, tuple], time_axis : int = 0, boundary : int = 0,
                    variable_names = ['u',], derivs = None, method = 'ANN', method_kwargs : dict = {},
                    max_deriv_order = 1, additional_tokens = [], coordinate_tensors = None,
                    memory_for_cache = 5, prune_domain : bool = False,
                    pivotal_tensor_label = None, pruner = None, threshold : float = 1e-2,
                    division_fractions = 3, rectangular : bool = True, data_fun_pow : int = 1):
        assert (isinstance(derivs, list) and isinstance(derivs[0], np.ndarray)) or derivs is None
        if isinstance(data, np.ndarray):
            data = [data,]

        set_grids = coordinate_tensors is not None
        set_grids_among_tokens = coordinate_tensors is not None
        if derivs is None:
            if len(data) != len(variable_names):
                print(len(data), len(variable_names))
                raise ValueError('Mismatching lengths of data tensors and the names of the variables')
        else:
            if not (len(data) == len(variable_names) == len(derivs)):
                raise ValueError('Mismatching lengths of data tensors, names of the variables and passed derivatives')
        data_tokens = []
        for data_elem_idx, data_tensor in enumerate(data):
            assert isinstance(data_tensor, np.ndarray), 'Input data must be in format of numpy ndarrays or iterable (list or tuple) of numpy arrays'
            entry = Input_data_entry(var_name = variable_names[data_elem_idx],
                                     data_tensor = data_tensor,
                                     coord_tensors = coordinate_tensors)
            derivs_tensor = derivs[data_elem_idx] if derivs is not None else None
            entry.set_derivatives(deriv_tensors = derivs_tensor, method = method, max_order = max_deriv_order,
                                  method_kwargs=method_kwargs)
            print(f'set grids parameter is {set_grids}')
            entry.use_global_cache(grids_as_tokens = set_grids_among_tokens,
                                   set_grids=set_grids, memory_for_cache=memory_for_cache, boundary=boundary)
            set_grids = False; set_grids_among_tokens = False

            entry_token_family = Token_family(entry.var_name, family_of_derivs = True)
            entry_token_family.set_status(unique_specific_token=False, unique_token_type=False,
                                 s_and_d_merged = False, meaningful = True,
                                 unique_for_right_part = True)
            entry_token_family.set_params(entry.names, OrderedDict([('power', (1, data_fun_pow))]),
                                          {'power' : 0}, entry.d_orders)
            entry_token_family.set_evaluator(simple_function_evaluator, [])

            print(entry_token_family.tokens)
            data_tokens.append(entry_token_family)
        if prune_domain:
            self.set_domain_pruning(pivotal_tensor_label, pruner, threshold, division_fractions, rectangular)

        if isinstance(additional_tokens, list):
            if not all([isinstance(tf, (Token_family, Prepared_tokens)) for tf in additional_tokens]):
                raise TypeError(f'Incorrect type of additional tokens: expected list or Token_family/Prepared_tokens - obj, instead got list of {type(additional_tokens[0])}')
        elif isinstance(additional_tokens, (Token_family, Prepared_tokens)):
            additional_tokens = [additional_tokens,]
        else:
            print(isinstance(additional_tokens, Prepared_tokens))
#            print(isinstance(additional_tokens, Custom_tokens))
#            print(issubclass(Custom_tokens, Prepared_tokens) or issubclass(Prepared_tokens, Custom_tokens))
            raise TypeError(f'Incorrect type of additional tokens: expected list or Token_family/Prepared_tokens - obj, instead got {type(additional_tokens)}')
        self.pool = TF_Pool(data_tokens + [tf if isinstance(tf, Token_family) else tf.token_family
                                      for tf in additional_tokens])
        print(f'The cardinality of defined token pool is {self.pool.families_cardinality()}')

    def fit(self, data : Union[np.ndarray, list, tuple], time_axis : int = 0, boundary : int = 0,
            equation_terms_max_number = 6, equation_factors_max_number = 1, variable_names = ['u',],
            eq_sparsity_interval = (1e-4, 2.5), derivs = None, max_deriv_order = 1,
            deriv_method = 'ANN', deriv_method_kwargs : dict = {},
            additional_tokens = [], coordinate_tensors = None, memory_for_cache = 5,
            prune_domain : bool = False, pivotal_tensor_label = None, pruner = None,
            threshold : float = 1e-2, division_fractions = 3, rectangular : bool = True,
            data_fun_pow : int = 1):
        '''

        Fit epde search algorithm to obtain differential equations, describing passed data.

        Parameters:
        -------------

        data : np.ndarray, list or tuple,
            Values of modeled variables. If the variable is single (i.e. deriving a single equation),
            it can be passed as the numpy.ndarray or as the list/tuple with a single element;
            multiple variables are not supported yet, use older interfaces.

        time_axis : int, optional,
            Index of the axis in the data, representing time. Default value: 0.

        boundary : int, optional,
            The number of grid nodes, considered as the boundary of the studied domain,
            that are excluded during equation search to avoid wrong results, linked to
            errors in derivatives estimation. Default value: 0.

        equation_terms_max_number : int, optional,
            The maximum number of terms, present in the derived equations. Default value: 6.

        equation_factors_max_number : int, optional,
            The maximum number of factors (token functions; real-valued coefficients are not counted here),
            present in terms of the equaton. Default value: 1.

        variable_names : list of str, optional,
            Names of the independent variables, passed into search mechanism. Length of the list must correspond
            to the number of np.ndarrays, sent with in ``data`` parameter. Defalut value: ``['u',]``, representing
            a single variable *u*.

        eq_sparsity_interval : tuple of two float values, optional,
            The left and right boundaries of interval with sparse regression values. Influence the number of
            terms in the equation. Default value: ``(1e-4, 2.5)``.

        derivs : list of lists of np.ndarrays or None, optional,
            Pre-computed values of derivatives. If ``None`` is passed, the derivatives are calculated in the
            method. Recommended to use, if the computations of derivatives take too long. For further information
            about using data, prepared in advance, check ``epde.prep.derivatives.Preprocess_derivatives`` function.
            Default value: None.

        max_deriv_order : int or tuple/list, optional,
            Highest order of calculated derivatives. Default value: 1.

        additional_tokens : None or list of ``Token_family`` or ``Prepared_tokens`` objects, optional
            Additional tokens, that would be used to construct the equations among the main variables and their
            derivatives. Objects of this list must be of type ``epde.interface.token_family.Token_family`` or
            of ``epde.interface.prepared_tokens.Prepared_tokens`` subclasses types. Default value: None.

        coordinate_tensors : list of np.ndarrays, optional
            Values of the coordinates on the grid nodes with studied functions values. In case of 1D-problem,
            that will be ``numpy.array``, while the parameter for higher dimensionality problems can be set from
            ``numpy.meshgrid`` function. With None, the tensors will be created as ranges with step of 1 between
            nodes. Defalut value: None.

        field_smooth : bool, optional
            Parameter, if the input variable fields shall be smoothed to avoid the errors. If the data is
            assumed to be noiseless, shall be set to False, otherwise - True. Default value: False.

        memory_for_cache : int or float, optional
            Limit for the cache (in fraction of the memory) for precomputed tensor values to be stored:
            if int, will be considered as the percentage of the entire memory, and if float,
            then as a fraction of memory.

        data_fun_pow : int
            Maximum power of token,

        '''

        self.create_pool(data = data, time_axis=time_axis, boundary=boundary, variable_names=variable_names,
                         derivs=derivs, method=deriv_method, method_kwargs=deriv_method_kwargs,
                         max_deriv_order=max_deriv_order, additional_tokens=additional_tokens,
                         coordinate_tensors=coordinate_tensors, memory_for_cache=memory_for_cache,
                         prune_domain=prune_domain, pivotal_tensor_label=pivotal_tensor_label,
                         pruner=pruner, threshold=threshold, division_fractions=division_fractions,
                         rectangular=rectangular, data_fun_pow=data_fun_pow)

        pop_constructor = operators.systems_population_constructor(pool = self.pool, terms_number = equation_terms_max_number,
                                                               max_factors_in_term=equation_factors_max_number,
                                                               eq_search_evo=self.director.constructor.strategy,
                                                               sparcity_interval = eq_sparsity_interval)

        self.moeadd_params['pop_constructor'] = pop_constructor
        self.optimizer = moeadd_optimizer(**self.moeadd_params)
        evo_operator = operators.sys_search_evolutionary_operator(operators.mixing_xover,
                                                                  operators.gaussian_mutation)
        self.optimizer.set_evolutionary(operator=evo_operator)
        best_obj = np.concatenate((np.ones([1,]),
                                  np.zeros(shape=len([1 for token_family in self.pool.families if token_family.status['meaningful']]))))
        self.optimizer.pass_best_objectives(*best_obj)
        self.optimizer.optimize(**self.moeadd_optimization_params)

        print('The optimization has been conducted.')
        self.search_conducted = True

    @property
    def equations_pareto_frontier(self):
        if not self.search_conducted:
            raise AttributeError('Pareto set of the best equations has not been discovered. Use ``self.fit`` method.')
        return self.optimizer.pareto_levels.levels

    def equation_search_results(self, only_print : bool = True, level_num = 1):
        if only_print:
            for idx in range(min(level_num, len(self.equations_pareto_frontier))):
                print('\n')
                print(f'{idx}-th non-dominated level')
                print('\n')
                [print(f'{solution.text_form} , with objective function values of {solution.obj_fun} \n')
                for solution in self.equations_pareto_frontier[idx]]
        else:
            return self.optimizer.pareto_levels.levels[:level_num]

    def solver_forms(self):
        '''
        Method returns solver forms of the equations on 1-st non-dominated levels in a form of Python list.
        '''
        sf = []
        for system in self.equations_pareto_frontier[0]:
            if len(system.structure) > 1:
                raise NotImplementedError('Currently, only "systems", containing a single equations, can be passed to solvers.')
            sf.append(system.structure[0].solver_form())
        return sf

    @property
    def cache(self):
        if global_var.grid_cache is not None:
            return global_var.grid_cache, global_var.tensor_cache
        else:
            return None, global_var.tensor_cache

#    def test_equation(self, equation_text_form : str):
#
#

max_stars_repo_name: vnleonenko/EPDE

max_issues_repo_name: vnleonenko/EPDE

max_forks_repo_name: vnleonenko/EPDE


text:
import cv2
import numpy as np

# Create window for image display
CASCADE_FN = "haarcascade_frontalface_default.xml"

# The scale used for face recognition.
# It is important as the face recognition algorithm works better on small images
# Also helps with removing faces that are too far away
RESIZE_SCALE = 3
RECTANGE_COLOUR = (117, 30, 104)
BOX_COLOR = (255, 255, 255)
THICKNESS = 2

def getFaceCoordinates(image):
  """Uses openCV to detect the face present in the input image.

  Returns a list of length 4, with the two corners of the rectangle that define
  the position of the face: [x1, y1, x2, y2], where (x1, y1) and (x2, y2)
  are the defining corners of the rectangle.
  """

  cascade = cv2.CascadeClassifier(CASCADE_FN)
  img_copy = cv2.resize(image, (image.shape[1]/RESIZE_SCALE,
                                image.shape[0]/RESIZE_SCALE))
  gray = cv2.cvtColor(img_copy, cv2.COLOR_BGR2GRAY)
  gray = cv2.equalizeHist(gray)
  rects = cascade.detectMultiScale(gray, 1.2, 3)

  # If there is no face or if we have more than 2 faces return None
  # because we do not deal with that yet
  if len(rects) != 1:
    return None

  r = rects[0]
  corners = [r[0], r[1], r[0] + r[2], r[1] + r[3]]

  return map((lambda x: RESIZE_SCALE * x), corners)

# Turns an image into a rbg image by replicating the 2d data on each of the axis.
def to_rgb1(im):
    w, h = im.shape
    ret = np.empty((w, h, 3), dtype=np.uint8)
    ret[:, :, 0] = im
    ret[:, :, 1] = im
    ret[:, :, 2] = im
    return ret

def drawFace(image, faceCoordinates, emotion, emotion_to_text, emotion_to_image=None):
  """ Draws the face information (together with emotion information) on the input image. """

  # Draw the face detection rectangles.
  cv2.rectangle(np.asarray(image),
                (faceCoordinates[0], faceCoordinates[1]),
                (faceCoordinates[2], faceCoordinates[3]),
                RECTANGE_COLOUR,
                thickness=THICKNESS)

  # Display the emotion on the webcam stream.
  if emotion is not None:
    cv2.putText(image,
                # Get the text associated with this emotion, but
                # if we do not have one just display the integer.
                emotion_to_text.get(emotion, ""),
                (faceCoordinates[0], faceCoordinates[2]),
                cv2.FONT_HERSHEY_SIMPLEX,
                2,
                BOX_COLOR,
                thickness=2)

    # Add a nice smiley to show the classification
    if emotion_to_image:
      smallImage = emotion_to_image[emotion]
      smallImage = cv2.resize(smallImage, (faceCoordinates[0], faceCoordinates[1]))
      smallImage = to_rgb1(smallImage)
      if smallImage.shape[0] > image.shape[0] or smallImage.shape[1] > image.shape[0]:
        return
      image[0:0+smallImage.shape[0], 0:0 + smallImage.shape[1]] = smallImage

def cropFace(image, faceCoordinates):
  return image[faceCoordinates[1]: faceCoordinates[3],
                faceCoordinates[0]: faceCoordinates[2]]


text:
"""
    This module demonstrates different filters. This implementation of those filters use GPU as a main resource.
"""
# * imports
from PIL import Image
import pycuda.driver as cuda
import pycuda.autoinit
from pycuda.compiler import SourceModule
import numpy
from datetime import datetime
from Utils import save_image, open_image
import pandas as pd

# * samples
DIR_IMG_OUTPUT = "ImgOutput"
SMALL_IMG = "tim-swaan-eOpewngf68w-unsplash_small.jpg"
MEDIUM_IMG = "tim-swaan-eOpewngf68w-unsplash_medium.jpg"
BIG_IMG = "tim-swaan-eOpewngf68w-unsplash_big.jpg"
BIGGEST_IMG = "tim-swaan-eOpewngf68w-unsplash_biggest.jpg"

EXCEL_NAME =  "OUTPUT_GPU.xlsx"

# * number of tests
N = 3

def grayscale_filter(image : Image):
    """
    This function filters an image using grayscale filter.

    About grayscale filter:
        Grayscale is a range of monochromatic shades from black to white.
        Therefore, a grayscale image contains only shades of gray and no color.

    Args:
        image (Image): Image to filter.

    Returns:
        Image: Image after grayscale filtering.
        float: time of memory allocation.
        float: time of operations in kernel.
        float: time of getting data back from gpu.
        float: total time of filtering.
    """

    total_timestamp_1 = datetime.now()

    image_pixels = numpy.array(image)
    image_pixels = image_pixels.astype(numpy.float32) # * Change data type

    alloc_timestamp_1 = datetime.now()

    device_image_pixels = cuda.mem_alloc(image_pixels.nbytes)  # * Allocate memory on the device, nbytes - number of bytes
    cuda.memcpy_htod(device_image_pixels, image_pixels)        # * Transfer the data to the GPU

    alloc_timestamp_2 = datetime.now()

    kernel_timestamp_1 = datetime.now()

    # * Kernel declaration - kernel grid and block size - https://en.wikipedia.org/wiki/Thread_block_(CUDA_programming)
    # ! 1D-indexing
    BLOCK_SIZE = 1024
    block = (1024, 1, 1)    # * 1024 threads in block (block size: 1024x1x1)
    check_size = numpy.int32(image.size[0] * image.size[1]) # width * height of image
    grid = (int(image.size[0] * image.size[1] / BLOCK_SIZE) + 1, 1, 1) # * grid size (number_of_blocks, 1, 1)

    # * Kernel definition - write CUDA C code
    kernel = """

    __global__ void grayscaleFilter( float * inputImage, int check ){

        int subpixelIndex = (threadIdx.x ) + blockDim.x * blockIdx.x ; // subpixelIndex = index of subpixel (pointing at red SP) -> [... 'r, g, b','r, g, b' ...]
        if(subpixelIndex *3 < check*3)
        {
            int value = 0.299 * inputImage[subpixelIndex*3] + 0.587 * inputImage[subpixelIndex*3 + 1] + 0.114 * inputImage[subpixelIndex*3 + 2];
            inputImage[subpixelIndex*3] = value;      // red subpixel
            inputImage[subpixelIndex*3+1] = value;    // green subpixel
            inputImage[subpixelIndex*3+2] = value;    // blue subpixel
        }
    }
    """

    # * Compile and get kernel function
    mod = SourceModule(kernel)

    func = mod.get_function("grayscaleFilter")
    func(device_image_pixels, check_size, block = block, grid = grid)

    kernel_timestamp_2 = datetime.now()

    # * Get back data from gpu
    back_data_timestamp_1 = datetime.now()
    grayscale_image_pixels = numpy.empty_like(image_pixels)         # * Return a new array with the same shape and type as a given array
    cuda.memcpy_dtoh(grayscale_image_pixels, device_image_pixels)   # * Fetch the data back from the GPU and display it
    grayscale_image_pixels = (numpy.uint8(grayscale_image_pixels))  # * Change data type
    back_data_timestamp_2 = datetime.now()

    grayscale_image = Image.fromarray(grayscale_image_pixels, mode ="RGB") # * Convert pixels to image

    total_timestamp_2 = datetime.now()

    # * Calculate total time operations
    alloc_time = (alloc_timestamp_2 - alloc_timestamp_1).total_seconds()*1000
    kernel_time = (kernel_timestamp_2 - kernel_timestamp_1).total_seconds()*1000
    back_data_time = (back_data_timestamp_2 - back_data_timestamp_1).total_seconds()*1000
    total_time = (total_timestamp_2 - total_timestamp_1).total_seconds()*1000

    return grayscale_image, alloc_time, kernel_time, back_data_time, total_time

def color_filter(image : Image, color_number: int):
    """
    This function filters an image using color filter.

    About color filter:
        This filter allows through only choosen (red/green/blue) color.

    Args:
        image (Image): Image to filter.
        color (int): Color of filter (0 - red, 1 - green, 2 - blue).

    Returns:
        Image: Image after choosen color filtering.
        float: time of memory allocation.
        float: time of operations in kernel.
        float: time of getting data back from gpu.
        float: total time of filtering.
    """

    total_timestamp_1 = datetime.now()

    image_pixels = numpy.array(image)
    image_pixels = image_pixels.astype(numpy.float32) # * Change data type

    alloc_timestamp_1 = datetime.now()

    device_image_pixels = cuda.mem_alloc(image_pixels.nbytes)  # * Allocate memory on the device
    cuda.memcpy_htod(device_image_pixels, image_pixels)        # * Transfer the data to the GPU

    alloc_timestamp_2 = datetime.now()

    kernel_timestamp_1 = datetime.now()

    # * Kernel declaration - kernel grid and block size
    BLOCK_SIZE = 1024
    block = (1024, 1, 1)
    check_size = numpy.int32(image.size[0] * image.size[1])
    grid = (int(image.size[0] * image.size[1] / BLOCK_SIZE) + 1, 1, 1)

    # * Kernel definition - write CUDA C code
    kernel = """

    __global__ void colorFilter( float * inputImage, int check, int color){

        int subpixelIndex = (threadIdx.x ) + blockDim.x * blockIdx.x; // subpixelIndex = index of subpixel (pointing at red SP) -> [... 'r, g, b','r, g, b' ...]
        if(subpixelIndex *3 < check*3)
        {
            if(color == 0) // red
            {
                inputImage[subpixelIndex*3] = inputImage[subpixelIndex*3];
                inputImage[subpixelIndex*3+1] = 0;
                inputImage[subpixelIndex*3+2] = 0;
            }
            else if(color == 1) //green
            {
                inputImage[subpixelIndex*3] = 0;
                inputImage[subpixelIndex*3+1] = inputImage[subpixelIndex*3+1];
                inputImage[subpixelIndex*3+2] = 0;
            }
            else if(color == 2) //blue
            {
                inputImage[subpixelIndex*3] = 0;
                inputImage[subpixelIndex*3+1] = 0;
                inputImage[subpixelIndex*3+2] = inputImage[subpixelIndex*3+2];
            }
        }
    }
    """

    color = color_number
    color = numpy.int32(color)

    # * Compile and get kernel function
    mod = SourceModule(kernel)
    func = mod.get_function("colorFilter")
    func(device_image_pixels, check_size, block = block, grid = grid)

    kernel_timestamp_2 = datetime.now()

    back_data_timestamp_1 = datetime.now()

    # * Get back data from gpu
    color_image_pixels = numpy.empty_like(image_pixels)         # * Return a new array with the same shape and type as a given array
    cuda.memcpy_dtoh(color_image_pixels, device_image_pixels)   # * Fetch the data back from the GPU and display it
    color_image_pixels = (numpy.uint8(color_image_pixels))      # * Change data type

    back_data_timestamp_2 = datetime.now()

    color_image = Image.fromarray(color_image_pixels, mode ="RGB") # * Convert pixels to image

    total_timestamp_2 = datetime.now()

    # * Calculate total time operations
    alloc_time = (alloc_timestamp_2 - alloc_timestamp_1).total_seconds()*1000
    kernel_time = (kernel_timestamp_2 - kernel_timestamp_1).total_seconds()*1000
    back_data_time = (back_data_timestamp_2 - back_data_timestamp_1).total_seconds()*1000
    total_time = (total_timestamp_2 - total_timestamp_1).total_seconds()*1000

    return color_image, alloc_time, kernel_time, back_data_time, total_time

def negative_filter(image: Image):
    """
    This function filters an image using negative filter.

    About negative filter:
        A positive image is a normal image. A negative image is a total inversion,
        in which light areas appear dark and vice versa. A negative color image is
        additionally color-reversed, with red areas appearing cyan, greens
        appearing magenta, and blues appearing yellow, and vice versa.

    Args:
        image (Image): Image to filter.

    Returns:
        Image: Image after negative filtering.
        float: time of memory allocation.
        float: time of operations in kernel.
        float: time of getting data back from gpu.
        float: total time of filtering.
    """

    total_timestamp_1 = datetime.now()

    image_pixels = numpy.array(image)
    image_pixels = image_pixels.astype(numpy.float32) # * Change data type

    alloc_timestamp_1 = datetime.now()

    device_image_pixels = cuda.mem_alloc(image_pixels.nbytes)  # * Allocate memory on the device
    cuda.memcpy_htod(device_image_pixels, image_pixels)        # * Transfer the data to the GPU

    alloc_timestamp_2 = datetime.now()

    kernel_timestamp_1 = datetime.now()

    # * Kernel declaration - kernel grid and block size
    BLOCK_SIZE = 1024
    block = (1024, 1, 1)
    check_size = numpy.int32(image.size[0] * image.size[1])
    grid = (int(image.size[0] * image.size[1] / BLOCK_SIZE) + 1, 1, 1)

    # * Kernel definition - write CUDA C code
    kernel = """

    __global__ void negativeFilter( float *inputImage, int check ){

        int subpixelIndex = (threadIdx.x ) + blockDim.x * blockIdx.x ; // subpixelIndex = index of subpixel (pointing at red SP) -> [... 'r, g, b','r, g, b' ...]

        if(subpixelIndex *3 < check*3)
        {
            inputImage[subpixelIndex*3] = 255 - inputImage[subpixelIndex*3];
            inputImage[subpixelIndex*3+1] = 255 - inputImage[subpixelIndex*3+1];
            inputImage[subpixelIndex*3+2] = 255 - inputImage[subpixelIndex*3+2];
        }
    }
    """

    # * Compile and get kernel function
    mod = SourceModule(kernel)

    func = mod.get_function("negativeFilter")
    func(device_image_pixels, check_size, block = block, grid = grid)

    kernel_timestamp_2 = datetime.now()

    back_data_timestamp_1 = datetime.now()

    # * Get back data from gpu
    negative_image_pixels = numpy.empty_like(image_pixels)         # * Return a new array with the same shape and type as a given array
    cuda.memcpy_dtoh(negative_image_pixels, device_image_pixels)   # * Fetch the data back from the GPU and display it
    negative_image_pixels = (numpy.uint8(negative_image_pixels))   # * Change data type

    back_data_timestamp_2 = datetime.now()

    negative_image = Image.fromarray(negative_image_pixels, mode ="RGB")    # * Convert pixels to image

    total_timestamp_2 = datetime.now()

    # * Calculate total time of operations
    alloc_time = (alloc_timestamp_2 - alloc_timestamp_1).total_seconds()*1000
    kernel_time = (kernel_timestamp_2 - kernel_timestamp_1).total_seconds()*1000
    back_data_time = (back_data_timestamp_2 - back_data_timestamp_1).total_seconds()*1000
    total_time = (total_timestamp_2 - total_timestamp_1).total_seconds()*1000

    return negative_image, alloc_time, kernel_time, back_data_time, total_time

def add_new_row_to_report(report : dict, test_img_name : str, filter_name : str,
 alloc_time: float, kernel_time: float, back_data_time: float, total_time: float):
    """
    This function adds given data to the dictionary (report).

    Args:
        report (dict): report that later will be converted to the excel.
        test_img_name (str): name of the image that is being tested.
        filter_name (str): name of the filter that is being tested.
        alloc_time (float): time of memory allocation.
        kernel_time (float): time of operations in kernel.
        back_data_time (float): time of getting data back from gpu.
        total_time (float): total time of filtering.

	Returns:
		dict: Report after update.
    """

    report['img_name'].append(test_img_name)
    report['filter_name'].append(filter_name)
    report['alloc_time'].append(alloc_time)
    report['kernel_time'].append(kernel_time)
    report['back_data_time'].append(back_data_time)
    report['total_time'].append(total_time)
    return report

def test_filter(test_img_name : str, report : dict, filter_name : str, filter_option : int):
    """
	This function

	Args:
		test_img_name (str): name of the iamge that is being tested.
        report (dict): report that later will be converted to the excel.
        filter_name (str): name of the filter that is being tested.
        filter_option (int): number that indicates the selected color filtering option.

	Returns:
		dict: Report after update.
	"""

    test_img = open_image(test_img_name)

    for i in range(0, N):
        if filter_name == 'grayscale':
            filtered_img_grayscale, alloc_time, kernel_time, back_data_time, total_time = grayscale_filter(test_img)
            report = add_new_row_to_report(report, test_img_name, filter_name, alloc_time, kernel_time, back_data_time, total_time)
            save_image(filtered_img_grayscale, f'./{DIR_IMG_OUTPUT}/grayscale_test_gpu.jpg')

        if filter_name == 'negative':
            filtered_img_negative, alloc_time, kernel_time, back_data_time, total_time = negative_filter(test_img)
            report = add_new_row_to_report(report, test_img_name, filter_name, alloc_time, kernel_time, back_data_time, total_time)
            save_image(filtered_img_negative, f'./{DIR_IMG_OUTPUT}/negative_test_gpu.jpg')

        if filter_name == 'color':
            filtered_img_color, alloc_time, kernel_time, back_data_time, total_time = color_filter(test_img, filter_option)
            report = add_new_row_to_report(report, test_img_name, filter_name, alloc_time, kernel_time, back_data_time, total_time)
            save_image(filtered_img_color, f'./{DIR_IMG_OUTPUT}/color_test_gpu.jpg')

    return report

if __name__ == "__main__":
    """
        SMALL_IMG - 640x427, 93.9 KB
        MEDIUM_IMG - 1920x1280, 838.8 KB
        BIG_IMG - 2400x1600, 1.24 MB
        BIGGEST_IMG - 5472x3648, 5.6 MB
    """

    report = {'img_name': [], 'filter_name': [], 'alloc_time': [], 'kernel_time': [], 'back_data_time': [], 'total_time': []}

    # * SMALL
    report = test_filter(test_img_name = SMALL_IMG, report = report, filter_name = 'grayscale', filter_option = None)
    report = test_filter(test_img_name = SMALL_IMG, report = report, filter_name = 'negative', filter_option = None)
    report = test_filter(test_img_name = SMALL_IMG, report = report, filter_name = 'color', filter_option = 0)
    # * MEDIUM
    report = test_filter(test_img_name = MEDIUM_IMG, report = report, filter_name = 'grayscale', filter_option = None)
    report = test_filter(test_img_name = MEDIUM_IMG, report = report, filter_name = 'negative', filter_option = None)
    report = test_filter(test_img_name = MEDIUM_IMG, report = report, filter_name = 'color', filter_option = 0)
    # * BIG
    report = test_filter(test_img_name = BIG_IMG, report = report, filter_name = 'grayscale', filter_option = None)
    report = test_filter(test_img_name = BIG_IMG, report = report, filter_name = 'negative', filter_option = None)
    report = test_filter(test_img_name = BIG_IMG, report = report, filter_name = 'color', filter_option = 0)
    # * BIGGEST
    report = test_filter(test_img_name = BIGGEST_IMG, report = report, filter_name = 'grayscale', filter_option = None)
    report = test_filter(test_img_name = BIGGEST_IMG, report = report, filter_name = 'negative', filter_option = None)
    report = test_filter(test_img_name = BIGGEST_IMG, report = report, filter_name = 'color', filter_option = 0)

    data = pd.DataFrame(report)
    data.to_excel(EXCEL_NAME)

max_stars_repo_name: Workata/Accelerated-Computing

max_issues_repo_name: Workata/Accelerated-Computing

max_forks_repo_name: Workata/Accelerated-Computing

reason: import numpy,import pycuda,from pycuda


text:
"""
    This file contains evaluation methods that take in a set of predicted labels
        and a set of ground truth labels and calculate precision, recall, accuracy, f1, and metrics @k
"""
from collections import defaultdict
import csv
import json
import numpy as np
import os
import sys

from sklearn.metrics import roc_curve, auc
from tqdm import tqdm

from constants import *
import datasets

def all_metrics(yhat, y, k=8, yhat_raw=None, calc_auc=True):
    """
        Inputs:
            yhat: binary predictions matrix
            y: binary ground truth matrix
            k: for @k metrics
            yhat_raw: prediction scores matrix (floats)
        Outputs:
            dict holding relevant metrics
    """
    names = ["acc", "prec", "rec", "f1"]

    #macro
    macro = all_macro(yhat, y)

    #micro
    ymic = y.ravel()
    yhatmic = yhat.ravel()
    micro = all_micro(yhatmic, ymic)

    metrics = {names[i] + "_macro": macro[i] for i in range(len(macro))}
    metrics.update({names[i] + "_micro": micro[i] for i in range(len(micro))})

    #AUC and @k
    if yhat_raw is not None and calc_auc:
        #allow k to be passed as int or list
        if type(k) != list:
            k = [k]
        for k_i in k:
            rec_at_k = recall_at_k(yhat_raw, y, k_i)
            metrics['rec_at_%d' % k_i] = rec_at_k
            prec_at_k = precision_at_k(yhat_raw, y, k_i)
            metrics['prec_at_%d' % k_i] = prec_at_k
            metrics['f1_at_%d' % k_i] = 2*(prec_at_k*rec_at_k)/(prec_at_k+rec_at_k)

        roc_auc = auc_metrics(yhat_raw, y, ymic)
        metrics.update(roc_auc)

    return metrics

def all_macro(yhat, y):
    return macro_accuracy(yhat, y), macro_precision(yhat, y), macro_recall(yhat, y), macro_f1(yhat, y)

def all_micro(yhatmic, ymic):
    return micro_accuracy(yhatmic, ymic), micro_precision(yhatmic, ymic), micro_recall(yhatmic, ymic), micro_f1(yhatmic, ymic)

#########################################################################
#MACRO METRICS: calculate metric for each label and average across labels
#########################################################################

def macro_accuracy(yhat, y):
    num = intersect_size(yhat, y, 0) / (union_size(yhat, y, 0) + 1e-10)
    return np.mean(num)

def macro_precision(yhat, y):
    num = intersect_size(yhat, y, 0) / (yhat.sum(axis=0) + 1e-10)
    return np.mean(num)

def macro_recall(yhat, y):
    num = intersect_size(yhat, y, 0) / (y.sum(axis=0) + 1e-10)
    return np.mean(num)

def macro_f1(yhat, y):
    prec = macro_precision(yhat, y)
    rec = macro_recall(yhat, y)
    if prec + rec == 0:
        f1 = 0.
    else:
        f1 = 2*(prec*rec)/(prec+rec)
    return f1

###################
# INSTANCE-AVERAGED
###################

def inst_precision(yhat, y):
    num = intersect_size(yhat, y, 1) / yhat.sum(axis=1)
    #correct for divide-by-zeros
    num[np.isnan(num)] = 0.
    return np.mean(num)

def inst_recall(yhat, y):
    num = intersect_size(yhat, y, 1) / y.sum(axis=1)
    #correct for divide-by-zeros
    num[np.isnan(num)] = 0.
    return np.mean(num)

def inst_f1(yhat, y):
    prec = inst_precision(yhat, y)
    rec = inst_recall(yhat, y)
    f1 = 2*(prec*rec)/(prec+rec)
    return f1

##############
# AT-K
##############

def recall_at_k(yhat_raw, y, k):
    #num true labels in top k predictions / num true labels
    sortd = np.argsort(yhat_raw)[:,::-1]
    topk = sortd[:,:k]

    #get recall at k for each example
    vals = []
    for i, tk in enumerate(topk):
        num_true_in_top_k = y[i,tk].sum()
        denom = y[i,:].sum()
        vals.append(num_true_in_top_k / float(denom))

    vals = np.array(vals)
    vals[np.isnan(vals)] = 0.

    return np.mean(vals)

def precision_at_k(yhat_raw, y, k):
    #num true labels in top k predictions / k
    sortd = np.argsort(yhat_raw)[:,::-1]
    topk = sortd[:,:k]

    #get precision at k for each example
    vals = []
    for i, tk in enumerate(topk):
        if len(tk) > 0:
            num_true_in_top_k = y[i,tk].sum()
            denom = len(tk)
            vals.append(num_true_in_top_k / float(denom))

    return np.mean(vals)

##########################################################################
#MICRO METRICS: treat every prediction as an individual binary prediction
##########################################################################

def micro_accuracy(yhatmic, ymic):
    return intersect_size(yhatmic, ymic, 0) / union_size(yhatmic, ymic, 0)

def micro_precision(yhatmic, ymic):
    return intersect_size(yhatmic, ymic, 0) / yhatmic.sum(axis=0)

def micro_recall(yhatmic, ymic):
    return intersect_size(yhatmic, ymic, 0) / ymic.sum(axis=0)

def micro_f1(yhatmic, ymic):
    prec = micro_precision(yhatmic, ymic)
    rec = micro_recall(yhatmic, ymic)
    if prec + rec == 0:
        f1 = 0.
    else:
        f1 = 2*(prec*rec)/(prec+rec)
    return f1

def auc_metrics(yhat_raw, y, ymic):
    if yhat_raw.shape[0] <= 1:
        return
    fpr = {}
    tpr = {}
    roc_auc = {}
    #get AUC for each label individually
    relevant_labels = []
    auc_labels = {}
    for i in range(y.shape[1]):
        #only if there are true positives for this label
        if y[:,i].sum() > 0:
            fpr[i], tpr[i], _ = roc_curve(y[:,i], yhat_raw[:,i])
            if len(fpr[i]) > 1 and len(tpr[i]) > 1:
                auc_score = auc(fpr[i], tpr[i])
                if not np.isnan(auc_score):
                    auc_labels["auc_%d" % i] = auc_score
                    relevant_labels.append(i)

    #macro-AUC: just average the auc scores
    aucs = []
    for i in relevant_labels:
        aucs.append(auc_labels['auc_%d' % i])
    roc_auc['auc_macro'] = np.mean(aucs)

    #micro-AUC: just look at each individual prediction
    yhatmic = yhat_raw.ravel()
    fpr["micro"], tpr["micro"], _ = roc_curve(ymic, yhatmic)
    roc_auc["auc_micro"] = auc(fpr["micro"], tpr["micro"])

    return roc_auc

########################
# METRICS BY CODE TYPE
########################

def results_by_type(Y, mdir, version='mimic3'):
    d2ind = {}
    p2ind = {}

    #get predictions for diagnoses and procedures
    diag_preds = defaultdict(lambda: set([]))
    proc_preds = defaultdict(lambda: set([]))
    preds = defaultdict(lambda: set())
    with open('%s/preds_test.psv' % mdir, 'r') as f:
        r = csv.reader(f, delimiter='|')
        for row in r:
            if len(row) > 1:
                for code in row[1:]:
                    preds[row[0]].add(code)
                    if code != '':
                        try:
                            pos = code.index('.')
                            if pos == 3 or (code[0] == 'E' and pos == 4):
                                if code not in d2ind:
                                    d2ind[code] = len(d2ind)
                                diag_preds[row[0]].add(code)
                            elif pos == 2:
                                if code not in p2ind:
                                    p2ind[code] = len(p2ind)
                                proc_preds[row[0]].add(code)
                        except:
                            if len(code) == 3 or (code[0] == 'E' and len(code) == 4):
                                if code not in d2ind:
                                    d2ind[code] = len(d2ind)
                                diag_preds[row[0]].add(code)
    #get ground truth for diagnoses and procedures
    diag_golds = defaultdict(lambda: set([]))
    proc_golds = defaultdict(lambda: set([]))
    golds = defaultdict(lambda: set())
    test_file = '%s/test_%s.csv' % (MIMIC_3_DIR, str(Y)) if version == 'mimic3' else '%s/test.csv' % MIMIC_2_DIR
    with open(test_file, 'r') as f:
        r = csv.reader(f)
        #header
        next(r)
        for row in r:
            codes = set([c for c in row[3].split(';')])
            for code in codes:
                golds[row[1]].add(code)
                try:
                    pos = code.index('.')
                    if pos == 3:
                        if code not in d2ind:
                            d2ind[code] = len(d2ind)
                        diag_golds[row[1]].add(code)
                    elif pos == 2:
                        if code not in p2ind:
                            p2ind[code] = len(p2ind)
                        proc_golds[row[1]].add(code)
                except:
                    if len(code) == 3 or (code[0] == 'E' and len(code) == 4):
                        if code not in d2ind:
                            d2ind[code] = len(d2ind)
                        diag_golds[row[1]].add(code)

    hadm_ids = sorted(set(diag_golds.keys()).intersection(set(diag_preds.keys())))

    ind2d = {i:d for d,i in d2ind.items()}
    ind2p = {i:p for p,i in p2ind.items()}
    type_dicts = (ind2d, ind2p)
    return diag_preds, diag_golds, proc_preds, proc_golds, golds, preds, hadm_ids, type_dicts

def diag_f1(diag_preds, diag_golds, ind2d, hadm_ids):
    num_labels = len(ind2d)
    yhat_diag = np.zeros((len(hadm_ids), num_labels))
    y_diag = np.zeros((len(hadm_ids), num_labels))
    for i,hadm_id in tqdm(enumerate(hadm_ids)):
        yhat_diag_inds = [1 if ind2d[j] in diag_preds[hadm_id] else 0 for j in range(num_labels)]
        gold_diag_inds = [1 if ind2d[j] in diag_golds[hadm_id] else 0 for j in range(num_labels)]
        yhat_diag[i] = yhat_diag_inds
        y_diag[i] = gold_diag_inds
    return micro_f1(yhat_diag.ravel(), y_diag.ravel())

def proc_f1(proc_preds, proc_golds, ind2p, hadm_ids):
    num_labels = len(ind2p)
    yhat_proc = np.zeros((len(hadm_ids), num_labels))
    y_proc = np.zeros((len(hadm_ids), num_labels))
    for i,hadm_id in tqdm(enumerate(hadm_ids)):
        yhat_proc_inds = [1 if ind2p[j] in proc_preds[hadm_id] else 0 for j in range(num_labels)]
        gold_proc_inds = [1 if ind2p[j] in proc_golds[hadm_id] else 0 for j in range(num_labels)]
        yhat_proc[i] = yhat_proc_inds
        y_proc[i] = gold_proc_inds
    return micro_f1(yhat_proc.ravel(), y_proc.ravel())

def metrics_from_dicts(preds, golds, mdir, ind2c):
    with open('%s/pred_100_scores_test.json' % mdir, 'r') as f:
        scors = json.load(f)

    hadm_ids = sorted(set(golds.keys()).intersection(set(preds.keys())))
    num_labels = len(ind2c)
    yhat = np.zeros((len(hadm_ids), num_labels))
    yhat_raw = np.zeros((len(hadm_ids), num_labels))
    y = np.zeros((len(hadm_ids), num_labels))
    for i,hadm_id in tqdm(enumerate(hadm_ids)):
        yhat_inds = [1 if ind2c[j] in preds[hadm_id] else 0 for j in range(num_labels)]
        yhat_raw_inds = [scors[hadm_id][ind2c[j]] if ind2c[j] in scors[hadm_id] else 0 for j in range(num_labels)]
        gold_inds = [1 if ind2c[j] in golds[hadm_id] else 0 for j in range(num_labels)]
        yhat[i] = yhat_inds
        yhat_raw[i] = yhat_raw_inds
        y[i] = gold_inds
    return yhat, yhat_raw, y, all_metrics(yhat, y, yhat_raw=yhat_raw, calc_auc=False)

def union_size(yhat, y, axis):
    #axis=0 for label-level union (macro). axis=1 for instance-level
    return np.logical_or(yhat, y).sum(axis=axis).astype(float)

def intersect_size(yhat, y, axis):
    #axis=0 for label-level union (macro). axis=1 for instance-level
    return np.logical_and(yhat, y).sum(axis=axis).astype(float)

def print_metrics(metrics):
    print()
    if "auc_macro" in metrics.keys():
        print("[MACRO] accuracy, precision, recall, f-measure, AUC")
        print("%.4f, %.4f, %.4f, %.4f, %.4f" % (metrics["acc_macro"], metrics["prec_macro"], metrics["rec_macro"], metrics["f1_macro"], metrics["auc_macro"]))
    else:
        print("[MACRO] accuracy, precision, recall, f-measure")
        print("%.4f, %.4f, %.4f, %.4f" % (metrics["acc_macro"], metrics["prec_macro"], metrics["rec_macro"], metrics["f1_macro"]))

    if "auc_micro" in metrics.keys():
        print("[MICRO] accuracy, precision, recall, f-measure, AUC")
        print("%.4f, %.4f, %.4f, %.4f, %.4f" % (metrics["acc_micro"], metrics["prec_micro"], metrics["rec_micro"], metrics["f1_micro"], metrics["auc_micro"]))
    else:
        print("[MICRO] accuracy, precision, recall, f-measure")
        print("%.4f, %.4f, %.4f, %.4f" % (metrics["acc_micro"], metrics["prec_micro"], metrics["rec_micro"], metrics["f1_micro"]))
    for metric, val in metrics.items():
        if metric.find("rec_at") != -1:
            print("%s: %.4f" % (metric, val))
    print()

if __name__ == "__main__":
    print(len(sys.argv))
    if len(sys.argv) < 4:
        print("usage: python " + str(os.path.basename(__file__) + " [train_dataset] [|Y| (as string)] [version (mimic2 or mimic3)] [model_dir （optional, default as current path)]")) #-HD
        sys.exit(0)
    elif len(sys.argv) == 4:
        sys.argv.append(os.getcwd())
    train_path, Y, version, mdir = sys.argv[1], sys.argv[2], sys.argv[3], sys.argv[4]
    ind2c, _ = datasets.load_full_codes(train_path, version=version)

    diag_preds, diag_golds, proc_preds, proc_golds, golds, preds, hadm_ids, type_dicts = results_by_type(Y, mdir, version)
    yhat, yhat_raw, y, metrics = metrics_from_dicts(preds, golds, mdir, ind2c)
    print_metrics(metrics)

    k = [5] if Y == '50' else [8,15]
    prec_at_8 = precision_at_k(yhat_raw, y, k=8)
    print("PRECISION@8: %.4f" % prec_at_8)
    prec_at_15 = precision_at_k(yhat_raw, y, k=15)
    print("PRECISION@15: %.4f" % prec_at_15)

    f1_diag = diag_f1(diag_preds, diag_golds, type_dicts[0], hadm_ids)
    f1_proc = proc_f1(proc_preds, proc_golds, type_dicts[1], hadm_ids)
    print("[BY CODE TYPE] f1-diag f1-proc")
    print("%.4f %.4f" % (f1_diag, f1_proc))

max_stars_repo_name: acadTags/caml-mimic

max_issues_repo_name: acadTags/caml-mimic

max_forks_repo_name: acadTags/caml-mimic


text:
import sys
import pandas as pd
import numpy as np

r_cols = ['user_id', 'song_link']
ratings = pd.read_csv('process/online.csv', sep='|', names=r_cols, usecols=range(2))

ratings = pd.DataFrame({'count' : ratings.groupby( [ "user_id", "song_link" ] ).size()}).reset_index()
ratings = ratings[ratings.user_id != sys.argv[2]]

if len(ratings[ratings.song_link == sys.argv[1]]) == 0:
	print '{"empty":"empty"}'

else :
	pivotTable = ratings.pivot_table(index=['user_id'], columns=['song_link'], values='count')

	pivotTable = pivotTable[np.isfinite(pivotTable[sys.argv[1]])]
	pivotTable = pivotTable.drop(sys.argv[1], 1).fillna(0)

	finalTable = pivotTable.mean(axis=0)
	finalTable = finalTable.sort_values(ascending=False)[:10]

	finalTable = finalTable.to_frame().reset_index()
	finalTable.columns = ["song_link", "ratings"]

	i_cols = ['user_id', 'song_link', 'name', 'artist', 'img']
	complete = pd.read_csv('process/online.csv', sep='|', names=i_cols, usecols=range(5))
	complete = complete.drop_duplicates(subset="song_link", keep="first")

	complete = complete.drop("user_id", 1)

	returnTable = pd.merge(finalTable, complete)
	returnTable = returnTable[returnTable.ratings != 0]
	returnTable = returnTable.drop("ratings", 1)

	print returnTable.reset_index().to_json(orient='records')


text:
import os
import sys
sys.path.insert(0, 'datasets')

import numpy as np

import torch
from torch.utils.data import DataLoader
from torchvision import datasets, transforms

from util.cutout import Cutout
from cifar import CustomCIFAR10

def create_loaders(args, hyper=False, root_dir='data/'):
    normalize = transforms.Normalize(mean=[x / 255.0 for x in [125.3, 123.0, 113.9]],
                                     std=[x / 255.0 for x in [63.0, 62.1, 66.7]])

    train_transform = transforms.Compose([])

    if not hyper:
        train_transform.transforms.append(transforms.ColorJitter(hue=args.hue,
            contrast=args.contrast, saturation=args.saturation,
            brightness=args.brightness))

    train_transform.transforms.append(transforms.ToTensor())
    train_transform.transforms.append(normalize)

    if not hyper:
        train_transform.transforms.append(Cutout(n_holes=args.cutholes, length=args.cutlength))

    test_transform = transforms.Compose([transforms.ToTensor(), normalize])

    if hyper:
        # Train set
        trainset = CustomCIFAR10(root=root_dir, train=True, download=True, transform=train_transform)
        num_train = int(np.floor((1-args.percent_valid) * len(trainset)))

        trainset.train_data = trainset.train_data[:num_train, :, :, :]
        trainset.train_labels = trainset.train_labels[:num_train]

        # Validation set
        valset = CustomCIFAR10(root=root_dir, train=True, download=True, transform=train_transform)
        valset.train_data = valset.train_data[num_train:, :, :, :]
        valset.train_labels = valset.train_labels[num_train:]

        # Test set
        testset = CustomCIFAR10(root=root_dir, train=False, download=True, transform=test_transform)

    else:
        # Train set
        trainset = datasets.CIFAR10(root=root_dir, train=True, download=True, transform=train_transform)
        num_train = int(np.floor((1-args.percent_valid) * len(trainset)))

        trainset.train_data = trainset.train_data[:num_train, :, :, :]
        trainset.train_labels = trainset.train_labels[:num_train]
        # Validation set
        valset = datasets.CIFAR10(root=root_dir, train=True, download=True, transform=train_transform)
        valset.train_data = valset.train_data[num_train:, :, :, :]
        valset.train_labels = valset.train_labels[num_train:]
        # Test set
        testset = datasets.CIFAR10(root=root_dir, train=False, download=True, transform=test_transform)

    train_loader = DataLoader(dataset=trainset, batch_size=args.batch_size, shuffle=True, pin_memory=False)
    valid_loader = DataLoader(dataset=valset, batch_size=args.batch_size, shuffle=True, pin_memory=False)
    test_loader = DataLoader(dataset=testset, batch_size=args.batch_size, pin_memory=False)

    return train_loader, valid_loader, test_loader

max_forks_repo_name: ThrunGroup/implicit-hyper-opt


text:
import numpy as np
import cv2

"""
List of functions:

* blank(shape, dtype = np.uint8, filler = '0')
    - creates a blank image (NumPy array) with either '0's or '1's

* simple_erode(img)
    - apply erosion with a 3x3 kernel and for 1 single iteration

* simple_dilate(img)
    - apply dilation with a 3x3 kernel and for 1 single iteration

* brightness_contrast(img, mult, add)
    - adjust the brightness and contrast of the image by multiplying and
      adding/subtracting values from the pixels;

* resize(img, width=None, height=None, inter=cv2.INTER_AREA)
    - resizes the image to the input height, width, or both

* getoutlines(img)
    - retrieve outlines (contours) of an input image

* order_points(pts)
    - arrange a list of corner points by coordinates, by the order of
      TOP-LEFT, TOP-RIGHT, BOTTOM-RIGHT, BOTTOM-LEFT

* perspective_transform(img, pts)
    - warps a deformed image into a straight, birds-eye view image,
      based on the corner points on the deformed image
"""

def blank(shape, dtype=np.uint8, filler="0"):
    """
    - creates a blank image (NumPy array) with either zeros or ones using
      built-in NumPy functions;
    - "shape" sets the shape of the blank image (array); takes tuples;
    - "dtype" sets the data type of the blank image, defaults to 8-bit
      unsigned integer; takes NumPy data-type objects;
    - "filler" sets the filler value of all pixels of the blank image, defaults
      to zeros; takes strings of either "0" or "1"
    """

    if filler == "0":
        blank = np.zeros(shape, dtype)

    elif filler == "1":
        blank = np.ones(shape, dtype)

    else:
        return "BAD FILLER VALUE; MUST BE STRINGS OF '0' OR '1'"

    return blank

def simple_erode(img):
    """
    - apply simple erosion to the input image using built-in OpenCV functions;
    - erosion kernel 3x3;
    - 1 single iteration
    """

    ekernel = np.ones((3, 3), np.uint8)
    eroded = cv2.erode(img, ekernel, iterations=1)

    return eroded

def simple_dilate(img):
    """
    - apply simple dilation to the input image using built-in OpenCV functions;
    - dilation kernel 3x3;
    - 1 single iteration
    """

    dkernel = np.ones((3, 3), np.uint8)
    dilated = cv2.dilate(img, dkernel, iterations=1)

    return dilated

def brightness_contrast(img, mult, add):
    """
    - allows brightness-contrast adjustment by multiplication and
      addition / subtraction of pixel values;
    - multiplication increases contrast while inevitably increasing overall
      brightness; set by "mult" parameter;
    - addition or subtraction increases or decreases value (brightness) of
      pixels; set by "add" parameter, use negative values for subtraction;
    """

    # multiply pixels by "mult", add by "add";
    # the multiplication will increase contrast, and adding/subtracting will
    # adjust the brightness
    adjusted = cv2.convertScaleAbs(img, alpha=float(mult), beta=float(add))

    return adjusted

# imutils resize() function.

"""
https://github.com/jrosebr1/imutils
MIT License.

Copyright (c) 2015-2016 Adrian Rosebrock, http://www.pyimagesearch.com

Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in
all copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
THE SOFTWARE.
"""

def resize(img, width=None, height=None, inter=cv2.INTER_AREA):
    """
    initialize the dimensions of the input image and obtain
    the image size
    """

    dim = None
    (h, w) = img.shape[:2]

    # if both the width and height are None, then return the
    # original image
    if width is None and height is None:
        return img

    # check to see if the width is None
    if width is None:
        # calculate the ratio of the height and construct the
        # dimensions
        r = height / float(h)
        dim = (int(w * r), height)

    # otherwise, the height is None
    else:
        # calculate the ratio of the width and construct the
        # dimensions
        r = width / float(w)
        dim = (width, int(h * r))

    # resize the image
    resized = cv2.resize(img, dim, interpolation=inter)

    # return the resized image
    return resized

def getoutlines(img):
    """
    - find the possible oulines of the input image using a built-in OpenCV
      function;
    - the "RETR_LIST" retreival mode returns a simple list of the found
      outlines;
    - the "cv2.CHAIN_APPROX_SIMPLE" approximation method returns coordinate
      points for the found outlines;
    - because the return of the contour function gives "contours", "heirarchy",
      we will only take the contours (outlines) for the current application
    """

    outlines = cv2.findContours(img, cv2.RETR_LIST, cv2.CHAIN_APPROX_SIMPLE)[0]
    return outlines

# 4-POINT PERSPECTIVE TRANSFORM FUNCTION SET.

"""
 - Call function "perspective_transform()";
 - Supply an input image, then 4 corner points;
 - The function will return the corrected image;
 - order_points() is a helper function
"""

def order_points(pts):
    """
    - returns a list of corner points in order
    - input "pts" will be a numpy array; will convert python list automatically;
    - "corners" list will go in the following order:
      0.    TOP-LEFT
      1.    TOP-RIGHT
      2. BOTTOM-RIGHT
      4. BOTTOM-LEFT
    """

    # automatically convert a python list to a numpy array
    if str(type(pts)) != "<class 'numpy.ndarray'>":
        pts = np.array(pts)

    # initialize an empty list to store the corner points
    corners = np.zeros((4, 2), dtype="float32")

    # figure out which set of coordinates are at which corner;

    # > the TOP-LEFT coordinates will have the smallest sum
    # > the BOTTOM-RIGHT coordinates will have the largest sum
    sums = pts.sum(axis=1)  # sum up all numbers horizontally
    corners[0] = pts[np.argmin(sums)]  # find out the TOP-LEFT coordinate
    corners[2] = pts[np.argmax(sums)]  # find out the BOTTOM-RIGHT coordinate

    # > the TOP-RIGHT coordinates will have the smallest difference
    # > the BOTTOM-LEFT coordinates will have the largest difference
    diffs = np.diff(pts, axis=1)
    corners[1] = pts[np.argmin(diffs)]
    corners[3] = pts[np.argmax(diffs)]

    return corners

def perspective_transform(img, pts):
    """
    - applies perspective transform to an image in order to straighten it,
      based on four given corner points
    - input "pts" will be a numpy array
    - returns a corrected image after applying perspective transform
    """

    # call the "order_points()" function to put corner points in order,
    # then assign each corner point to its respective variable
    corners_old = order_points(pts)
    tl, tr, br, bl = corners_old

    # calculate the WIDTH of the corrected image as follows:
    # > find the distance between the top points and the bottom points;
    # > find the maximum of the two distances, make it the new width.
    distT = np.sqrt(((tr[0] - tl[0]) ** 2) + ((tr[1] - tl[1]) ** 2))
    distB = np.sqrt(((br[0] - bl[0]) ** 2) + ((br[1] - bl[1]) ** 2))
    maxW = max(int(distT), int(distB))

    # calculate the HEIGHT of the corrected image as follows:
    # > find the distance between the left points and the right points;
    # > find the maximum of the two distances, make it the new height.
    distL = np.sqrt(((tl[0] - bl[0]) ** 2) + ((tl[1] - bl[1]) ** 2))
    distR = np.sqrt(((tr[0] - br[0]) ** 2) + ((tr[1] - br[1]) ** 2))
    maxH = max(int(distL), int(distR))

    # define corner points for the corrected image based on calculations
    # done above;
    # the same order of corners are followed
    corners_corrected = np.array(
        [[0, 0], [maxW - 1, 0], [maxW - 1, maxH - 1], [0, maxH - 1]], dtype="float32"
    )

    # using OpenCV, calculate a transform matrix then apply it to create
    # the corrected image
    matrix = cv2.getPerspectiveTransform(corners_old, corners_corrected)
    img_corrected = cv2.warpPerspective(img, matrix, (maxW, maxH))

    return img_corrected

max_stars_repo_name: hacksider/pyimgscan

max_issues_repo_name: hacksider/pyimgscan

max_forks_repo_name: hacksider/pyimgscan


text:
import TccLib

import numpy as np

# Import gmsh created mesh
mesh = TccLib.Mesh("Poisson_Dirichlet")

# Show mesh geometry
mesh.show_geometry(names=False, rainbow=True, save=True)

# Define parameters
T_0 = 0.
T_L = 0.
k = 5
Q = 40

# Define analytic comparison expression
analytic_expression = lambda x: Q/(2.*k) * (-x**2 + mesh.length_y*x) + (T_L-T_0)/mesh.length_y * x + T_0

# Define boundary conditions
boundary_conditions_values = {
    "north": T_L,
    "south": T_0,
}

boundary_conditions_types = {
    "north": True,
    "south": True,
}

# Get vectors of nodes values
xy_indices, xy_values, xy_types = TccLib.build_boundary_conditions(mesh, boundary_conditions_values,
                                                                   boundary_conditions_types)

# Set Boundary Conditions
mesh.new_boundary_condition("space", point_index=xy_indices, values=xy_values,
                            type_of_boundary=xy_types)
mesh.new_boundary_condition("time", point_index=xy_indices, values=xy_values,
                            type_of_boundary=xy_types)

# Solve for permanent solution
temperature_perm = TccLib.solve_poisson(mesh, permanent_solution=True, k_coef=k, q=Q)

# Show results in 3D graph
# mesh.show_3d_solution(temperature_perm, view_from_above=False, axis_labels={"z": "T(°C)"})

# Define x vector of positions
x_vector = np.linspace(min(mesh.y), max(mesh.y), 100)

# Find values of property in the mesh at a determined set position for every value in the vector of x
x_position = mesh.length_x/2.
y_vector = [mesh.get_interpolated_value([x_position, x], temperature_perm) for x in x_vector]

# Show comparison graph
TccLib.util.show_comparison(x_vector, analytic_expression, y_vector, numeric_label="Solução Numérica",
                            analytic_label="Solução Analítica", title="Equação de Poisson Permanente",
                            x_label="Posição no Eixo y (m)", y_label="Valor da Temperatura (°C)",
                            save_file_as="{0}_permanent_validation".format(mesh.name))

# Solve for transient solution
temperature_trans = TccLib.solve_poisson(mesh, permanent_solution=False, k_coef=k, stop_rule=1e-5, q=Q,
                                         return_history=True)

# Show results in 3D graph
mesh.show_animated_3d_solution(temperature_trans, dt=mesh.default_dt, axis_labels={"z": "T(°C)"})

# Get a small dictionary with each value state with timestamps as keys
small_dict = temperature_trans.reduced_dict_log(5)

# Find values of property in the mesh at a determined set position for every value in the vector of x at each timestamp
z_vectors = {key: [mesh.get_interpolated_value([x_position, x], small_dict[key]) for x in x_vector]
             for key in small_dict}

# Show comparison graph
TccLib.util.show_comparison(x_vector, analytic_expression, z_vectors, numeric_label="Solução Numérica",
                            analytic_label="Solução Analítica", title="Equação de Poisson Transiente",
                            x_label="Posição no Eixo y (m)", y_label="Valor da Temperatura (°C)",
                            save_file_as="{0}_transient_validation".format(mesh.name))

max_stars_repo_name: encarvlucas/LucasCarvalhoTCC

max_issues_repo_name: encarvlucas/LucasCarvalhoTCC

max_forks_repo_name: encarvlucas/LucasCarvalhoTCC


text:
from sympy import ( symbols, solve, diff, integrate, exp, sqrt, lambdify, pprint )
import matplotlib.pyplot as plt
import numpy as np

# The difference between the true value of an integral and the value given by the trapezoidal rule or Simpson's rule is known as the error.
# In numerical analysis, the error is studied to determine how large n must be for the error to be smaller than some specified amount.
#
# For both rules, the error is inversely proportional to a power of n, the number of subdivisions.

# In other words, the error is roughly StartFraction k / n^p
# where k is a constant that depends on the function and the interval, and p is a power that depends only on the method used.
# With a little experimentation, you can find out what the power p is for the trapezoidal rule and for Simpson's rule.

def simpsons_rule( f, a, b, n ):
	area = 0

	steps = np.linspace( a, b, n + 1, endpoint = True )

	for i in range( 0, n + 1 ):

		v = f.subs( { x: steps[ i ] } )

		if i == 0 or i == n:
			area += v
		elif i % 2 == 0:
			area += 2 * v
		else:
			area += 4 * v

	return area * ( b - a ) / ( 3*n )

# Function

x = symbols( 'x' )
F = 11*x**4

# What does the petal distribution look like?

g_xlim = [ -.01, 2 ]

lam_p = lambdify( x, F, np )
x_vals = np.linspace( g_xlim[0], g_xlim[1], 1000, endpoint=True )
y_vals = lam_p( x_vals )
plt.plot( x_vals, y_vals )

x_min, x_max = 0, 1

plt.vlines( x = x_min, ymin = 0, ymax = F.subs( { x: x_min } ), color = 'Black', zorder = 1 )
plt.vlines( x = x_max, ymin = 0, ymax = F.subs( { x: x_max } ), color = 'Red', zorder = 1 )

plt.show()

#  Find the exact value of
a, b = 0, 1
area = integrate( F, ( x, a, b ) ).evalf()
area

# Approximate the integral in part a using Simpson's rule with n = 4, 8, 16, and 32
n = 4
n4 = round( simpsons_rule( F, a, b, n ), 7 )
n4
e4 = round( abs( area - n4 ), 7 )
'{:.7f}'.format( abs( n4 - area ) )

n = 8
n8 = round( simpsons_rule( F, a, b, n ), 7 )
n8
e8 = abs( area - n8 )
'{:.7f}'.format( abs( n8 - area ) )

n = 16
n16 = round( simpsons_rule( F, a, b, n ), 7 )
n16
e16 = abs( n16 - area )
'{:.7f}'.format( e16 )

n = 32
n32 = round( simpsons_rule( F, a, b, n ), 7 )
n32
e32 =  abs( n32 - area )
'{:.7f}'.format( e32 )

# c
#  If the error is k/n**p, then the error times n Superscript p should be approximately a constant. Multiply the errors in part b times n**p for p=1, 2, etc.,
#  until you find a power p yielding approximately the same answer for all four values of n.

p = 4

e4 * 4 ** p
e8 * 8 ** p
e16 * 16 ** p
e32 * 32 ** p

max_stars_repo_name: bmoretz/Python-Playground

max_issues_repo_name: bmoretz/Python-Playground

max_forks_repo_name: bmoretz/Python-Playground

reason: import numpy,from sympy