text:
import os
import os.path
from pprint import pprint
import configparser
import csv
import fiona
import numpy as np
import random
import glob

from shapely.geometry import shape, Point, LineString, Polygon, MultiPolygon, mapping
from shapely.ops import unary_union, cascaded_union
from pyproj import Proj, transform
from rtree import index

from collections import OrderedDict, defaultdict

import osmnx as ox, networkx as nx

CONFIG = configparser.ConfigParser()
CONFIG.read(os.path.join(os.path.dirname(__file__), 'script_config.ini'))
BASE_PATH = CONFIG['file_locations']['base_path']

#####################################
# setup file locations and data files
#####################################

DATA_RAW = os.path.join(BASE_PATH, 'raw')
DATA_INTERMEDIATE = os.path.join(BASE_PATH, 'intermediate')

#####################################
# READ LOOK UP TABLE (LUT) DATA
#####################################

def read_pcd_to_exchange_lut():
    """
    Produces all unique postcode-to-exchange combinations from available data, including:

    'January 2013 PCP to Postcode File Part One.csv'
    'January 2013 PCP to Postcode File Part Two.csv'
    'pcp.to.pcd.dec.11.one.csv'
    'pcp.to.pcd.dec.11.two.csv'
    'from_tomasso_valletti.csv'

    Data Schema
    ----------
    * exchange_id: 'string'
        Unique Exchange ID
    * postcode: 'string'
        Unique Postcode

    Returns
    -------
    pcd_to_exchange_data: List of dicts
    """
    pcd_to_exchange_data = []

    with open(os.path.join(DATA_RAW, 'network_hierarchy_data', 'January 2013 PCP to Postcode File Part One.csv'), 'r', encoding='utf8', errors='replace') as system_file:
        reader = csv.reader(system_file)
        for skip in range(11):
            next(reader)
        for line in reader:
            pcd_to_exchange_data.append({
                'exchange_id': line[0],
                'postcode': line[1].replace(" ", "")
            })

    with open(os.path.join(DATA_RAW, 'network_hierarchy_data', 'January 2013 PCP to Postcode File Part One.csv'), 'r',  encoding='utf8', errors='replace') as system_file:
        reader = csv.reader(system_file)
        for skip in range(11):
            next(reader)
        for line in reader:
            pcd_to_exchange_data.append({
                'exchange_id': line[0],
                'postcode': line[1].replace(" ", "")
            })

    with open(os.path.join(DATA_RAW, 'network_hierarchy_data', 'pcp.to.pcd.dec.11.one.csv'), 'r',  encoding='utf8', errors='replace') as system_file:
        reader = csv.reader(system_file)
        for skip in range(11):
            next(reader)
        for line in reader:
            pcd_to_exchange_data.append({
                'exchange_id': line[0],
                'postcode': line[1].replace(" ", "")
            })

    with open(os.path.join(DATA_RAW, 'network_hierarchy_data', 'pcp.to.pcd.dec.11.two.csv'), 'r', encoding='utf8', errors='replace') as system_file:
        reader = csv.reader(system_file)
        next(reader)
        for line in reader:
            pcd_to_exchange_data.append({
                'exchange_id': line[0],
                'postcode': line[1].replace(" ", "")
            })

    with open(os.path.join(DATA_RAW, 'network_hierarchy_data', 'from_tomasso_valletti.csv'), 'r', encoding='utf8', errors='replace') as system_file:
        reader = csv.reader(system_file)
        next(reader)
        for line in reader:
            pcd_to_exchange_data.append({
                'exchange_id': line[0],
                'postcode': line[1].replace(" ", "")
            })

    ### find unique values in list of dicts
    return list({pcd['postcode']:pcd for pcd in pcd_to_exchange_data}.values())

def read_postcode_areas():

    """
    Reads all postcodes shapes, removing vertical postcodes, and merging with closest neighbour.

    Data Schema
    -----------
    * POSTCODE: 'string'
        Unique Postcode

    Returns
    -------
    postcode_areas = list of dicts
    """

    postcode_areas = []

    pathlist = glob.iglob(os.path.join(DATA_RAW, 'codepoint', 'codepoint-poly_2429451') + '/**/*.shp', recursive=True)

    for path in pathlist:

        # Initialze Rtree
        idx = index.Index()

        print(path)
        with fiona.open(path, 'r') as source:

            # Store shapes in Rtree
            for src_shape in source:
                idx.insert(int(src_shape['id']), shape(src_shape['geometry']).bounds, src_shape)

            # Split list in regular and vertical postcodes
            postcodes = {}
            vertical_postcodes = {}

            for x in source:

                x['properties']['POSTCODE'] = x['properties']['POSTCODE'].replace(" ", "")
                if x['properties']['POSTCODE'].startswith('V'):
                    vertical_postcodes[x['id']] = x
                else:
                    postcodes[x['id']] = x

            for key, f in vertical_postcodes.items():

                vpost_geom = shape(f['geometry'])
                best_neighbour = {'id': 0, 'intersection': 0}

                # Find best neighbour
                for n in idx.intersection((vpost_geom.bounds), objects=True):
                    if shape(n.object['geometry']).intersection(vpost_geom).length > best_neighbour['intersection'] and n.object['id'] != f['id']:
                        best_neighbour['id'] = n.object['id']
                        best_neighbour['intersection'] = shape(n.object['geometry']).intersection(vpost_geom).length

                # Merge with best neighbour
                neighbour = postcodes[best_neighbour['id']]
                merged_geom = unary_union([shape(neighbour['geometry']), vpost_geom])

                merged_postcode = {
                    'id': neighbour['id'].replace(" ", ""),
                    'properties': neighbour['properties'],
                    'geometry': mapping(merged_geom)
                }

                try:
                    postcodes[merged_postcode['id']] = merged_postcode
                except:
                    raise Exception

            for key, p in postcodes.items():
                p.pop('id')
                postcode_areas.append(p)

    return postcode_areas

def read_exchanges():

    """
    Reads in exchanges from 'final_exchange_pcds.csv'.

    Data Schema
    ----------
    * id: 'string'
        Unique Exchange ID
    * Name: 'string'
        Unique Exchange Name
    * pcd: 'string'
        Unique Postcode
    * Region: 'string'
        Region ID
    * County: 'string'
        County IS

    Returns
    -------
    exchanges: List of dicts
    """

    exchanges = []

    with open(os.path.join(DATA_RAW, 'layer_2_exchanges', 'final_exchange_pcds.csv'), 'r') as system_file:
        reader = csv.reader(system_file)
        next(reader)

        for line in reader:
            exchanges.append({
                'type': "Feature",
                'geometry': {
                    "type": "Point",
                    "coordinates": [float(line[5]), float(line[6])]
                },
                'properties': {
                    'id': 'exchange_' + line[1],
                    'Name': line[2],
                    'pcd': line[0],
                    'Region': line[3],
                    'County': line[4]
                }
            })

    return exchanges

#####################################
# PROCESS NETWORK HIERARCHY
#####################################

def add_exchange_id_to_postcode_areas(exchanges, postcode_areas, exchange_to_postcode):

    """
    Either uses known data or estimates which exchange each postcode is likely attached to.

    Arguments
    ---------

    * exchanges: 'list of dicts'
        List of Exchanges from read_exchanges()
    * postcode_areas: 'list of dicts'
        List of Postcode Areas from read_postcode_areas()
    * exchange_to_postcode: 'list of dicts'
        List of Postcode to Exchange data procudes from read_pcd_to_exchange_lut()

    Returns
    -------
    postcode_areas: 'list of dicts'
    """
    idx_exchanges = index.Index()
    lut_exchanges = {}

    # Read the exchange points
    for idx, exchange in enumerate(exchanges):

        # Add to Rtree and lookup table
        idx_exchanges.insert(idx, tuple(map(int, exchange['geometry']['coordinates'])) + tuple(map(int, exchange['geometry']['coordinates'])), exchange['properties']['id'])
        lut_exchanges[exchange['properties']['id']] = {
            'Name': exchange['properties']['Name'],
            'pcd': exchange['properties']['pcd'].replace(" ", ""),
            'Region': exchange['properties']['Region'],
            'County': exchange['properties']['County'],
        }

    # Read the postcode-to-cabinet-to-exchange lookup file
    lut_pcb2cab = {}

    for idx, row in enumerate(exchange_to_postcode):
        lut_pcb2cab[row['postcode']] = row['exchange_id']

    # Connect each postcode area to an exchange
    for postcode_area in postcode_areas:

        postcode = postcode_area['properties']['POSTCODE']

        if postcode in lut_pcb2cab:

            # Postcode-to-cabinet-to-exchange association
            postcode_area['properties']['EX_ID'] = 'exchange_' + lut_pcb2cab[postcode]
            postcode_area['properties']['EX_SRC'] = 'EXISTING POSTCODE DATA'

        else:

            # Find nearest exchange
            nearest = [n.object for n in idx_exchanges.nearest((shape(postcode_area['geometry']).bounds), 1, objects=True)]
            postcode_area['properties']['EX_ID'] = nearest[0]
            postcode_area['properties']['EX_SRC'] = 'ESTIMATED NEAREST'

        # Match the exchange ID with remaining exchange info
        if postcode_area['properties']['EX_ID'] in lut_exchanges:
            postcode_area['properties']['EX_NAME'] = lut_exchanges[postcode_area['properties']['EX_ID']]['Name']
            postcode_area['properties']['EX_PCD'] = lut_exchanges[postcode_area['properties']['EX_ID']]['pcd']
            postcode_area['properties']['EX_REGION'] = lut_exchanges[postcode_area['properties']['EX_ID']]['Region']
            postcode_area['properties']['EX_COUNTY'] = lut_exchanges[postcode_area['properties']['EX_ID']]['County']
        else:
            postcode_area['properties']['EX_NAME'] = ""
            postcode_area['properties']['EX_PCD'] = ""
            postcode_area['properties']['EX_REGION'] = ""
            postcode_area['properties']['EX_COUNTY'] = ""

    return postcode_areas

def generate_exchange_area(exchanges, merge=True):

    exchanges_by_group = defaultdict(list)

    # Loop through all exchanges
    print('generate_exchange_area - Group polygons by exchange ID')
    for f in exchanges:

        # Convert Multipolygons to list of polygons
        if (isinstance(shape(f['geometry']), MultiPolygon)):
            polygons = [p.buffer(0) for p in shape(f['geometry'])]
        else:
            polygons = [shape(f['geometry'])]

        exchanges_by_group[f['properties']['EX_ID']].extend(polygons)

    # Write Multipolygons per exchange
    print('generate_exchange_area - Generate multipolygons')
    exchange_areas = []
    for exchange, area in exchanges_by_group.items():

        exchange_multipolygon = MultiPolygon(area)
        exchange_areas.append({
            'type': "Feature",
            'geometry': mapping(exchange_multipolygon),
            'properties': {
                'id': exchange
            }
        })

    if merge:
        print('generate_exchange_area - Merge multipolygons into singlepolygons')
        # Merge MultiPolygons into single Polygon
        removed_islands = []
        for area in exchange_areas:

            # Avoid intersections
            geom = shape(area['geometry']).buffer(0)
            cascaded_geom = unary_union(geom)

            # Remove islands
            # Keep polygon with largest area
            # Add removed islands to a list so that they
            # can be merged in later
            if (isinstance(cascaded_geom, MultiPolygon)):
                for idx, p in enumerate(cascaded_geom):
                    if idx == 0:
                        geom = p
                    elif p.area > geom.area:
                        removed_islands.append(geom)
                        geom = p
                    else:
                        removed_islands.append(p)
            else:
                geom = cascaded_geom

            # Write exterior to file as polygon
            exterior = Polygon(list(geom.exterior.coords))

            # Write to output
            area['geometry'] = mapping(exterior)

        # Add islands that were removed because they were not
        # connected to the main polygon and were not recovered
        # because they were on the edge of the map or inbetween
        # exchanges :-). Merge to largest intersecting exchange area.
        print('generate_exchange_area - Process removed islands')
        idx_exchange_areas = index.Index()
        for idx, exchange_area in enumerate(exchange_areas):
            idx_exchange_areas.insert(idx, shape(exchange_area['geometry']).bounds, exchange_area)
        for island in removed_islands:
            intersections = [n for n in idx_exchange_areas.intersection((island.bounds), objects=True)]

            if len(intersections) > 0:
                for idx, intersection in enumerate(intersections):
                    if idx == 0:
                        merge_with = intersection
                    elif shape(intersection.object['geometry']).intersection(island).length > shape(merge_with.object['geometry']).intersection(island).length:
                        merge_with = intersection

                merged_geom = merge_with.object
                merged_geom['geometry'] = mapping(shape(merged_geom['geometry']).union(island))
                idx_exchange_areas.delete(merge_with.id, shape(merge_with.object['geometry']).bounds)
                idx_exchange_areas.insert(merge_with.id, shape(merged_geom['geometry']).bounds, merged_geom)

        exchange_areas = [n.object for n in idx_exchange_areas.intersection(idx_exchange_areas.bounds, objects=True)]

    return exchange_areas

def read_exchange_area():
    with fiona.open(os.path.join(DATA_INTERMEDIATE, '_exchange_areas.shp'), 'r') as source:
        return [exchange for exchange in source]

def write_shapefile(data, path):

    # Translate props to Fiona sink schema
    prop_schema = []
    for name, value in data[0]['properties'].items():
        fiona_prop_type = next((fiona_type for fiona_type, python_type in fiona.FIELD_TYPES_MAP.items() if python_type == type(value)), None)
        prop_schema.append((name, fiona_prop_type))

    sink_driver = 'ESRI Shapefile'
    sink_crs = {'init': 'epsg:27700'}
    sink_schema = {
        'geometry': data[0]['geometry']['type'],
        'properties': OrderedDict(prop_schema)
    }

    # Write all elements to output file
    with fiona.open(os.path.join(DATA_INTERMEDIATE, path), 'w', driver=sink_driver, crs=sink_crs, schema=sink_schema) as sink:
        for feature in data:
            sink.write(feature)

#####################################
# APPLY METHODS
#####################################

if __name__ == "__main__":

    SYSTEM_INPUT = os.path.join('data', 'digital_comms', 'raw')

    if not os.path.isfile(os.path.join(DATA_INTERMEDIATE, '_exchange_areas.shp')):

        # Read LUTs
        print('read_pcd_to_exchange_lut')
        lut_pcd_to_exchange = read_pcd_to_exchange_lut()

        print('read postcode_areas')
        geojson_postcode_areas = read_postcode_areas()

        # Write
        print('write postcode_areas')
        write_shapefile(geojson_postcode_areas, '_postcode_areas.shp')

        print('read exchanges')
        geojson_layer2_exchanges = read_exchanges()

        # Process/Estimate network hierarchy
        print('add exchange id to postcode areas')
        geojson_postcode_areas = add_exchange_id_to_postcode_areas(geojson_layer2_exchanges, geojson_postcode_areas, lut_pcd_to_exchange)

        print('generate exchange areas')
        geojson_exchange_areas = generate_exchange_area(geojson_postcode_areas)

        # Write
        print('write exchange_areas')
        write_shapefile(geojson_exchange_areas, '_exchange_areas.shp')

    exchange_areas = read_exchange_area()
    #[print(exchange['properties']['id']) for exchange in exchange_areas]

    print('exchange_EAARR')
    print('exchange_EABTM')
    print('exchange_EABWL')
    print('exchange_EACAM')
    print('exchange_EACFH')
    print('exchange_EACOM')
    print('exchange_EACRH')
    print('exchange_EACTM')
    print('exchange_EAESW')
    print('exchange_EAFUL')
    print('exchange_EAGIR')
    print('exchange_EAHIS')
    print('exchange_EAHST')
    print('exchange_EALNT')
    print('exchange_EAMAD')
    print('exchange_EAMBN')
    print('exchange_EASCI')
    print('exchange_EASIX')
    print('exchange_EASST')
    print('exchange_EASWV')
    print('exchange_EATEV')
    print('exchange_EATRU')
    print('exchange_EAWLM')
    print('exchange_EAWTB')

max_stars_repo_name: mchakli/test

max_issues_repo_name: mchakli/test

max_forks_repo_name: mchakli/test


text:
from __future__ import absolute_import, print_function, division

import unittest
import dill
import glob
import os

import numpy as np

from pysaliency.utils import LazyList, TemporaryDirectory, Cache, get_minimal_unique_filenames
from test_helpers import TestWithData

def test_minimal_unique_filenames():
    assert get_minimal_unique_filenames(['a/b/c.d']) == ['c.d']

    filenames = [
        'a/b/c/d.e',
        'a/b/c/f.g',
        'a/b/c/h.i',
    ]

    assert get_minimal_unique_filenames(filenames) == ['d.e', 'f.g', 'h.i']

    filenames.append('a/b/C/j.k')
    assert get_minimal_unique_filenames(filenames) == ['c/d.e', 'c/f.g', 'c/h.i', 'C/j.k']

class TestLazyList(TestWithData):
    def test_lazy_list(self):
        calls = []

        def gen(i):
            calls.append(i)
            print ('calling with {} yielding {}'.format(i, i**2))
            return i**2

        length = 20

        l = LazyList(gen, length)
        self.assertEqual(len(l), length)

        for i in range(length):
            self.assertEqual(l[i], i**2)

        self.assertEqual(calls, list(range(length)))

    def test_pickle_no_cache(self):
        def gen(i):
            print ('calling with {} yielding {}'.format(i, i**2))
            return i**2

        length = 20
        l = LazyList(gen, length)

        l = self.pickle_and_reload(l, pickler=dill)

        self.assertEqual(l._cache, {})
        self.assertEqual(list(l), [i**2 for i in range(length)])

    def test_pickle_with_cache(self):
        def gen(i):
            print ('calling with {} yielding {}'.format(i, i**2))
            return i**2

        length = 20
        l = LazyList(gen, length, pickle_cache=True)

        list(l)  # make sure all list items are generated

        l = self.pickle_and_reload(l, pickler=dill)

        self.assertEqual(l._cache, {i: i**2 for i in range(length)})
        self.assertEqual(list(l), [i**2 for i in range(length)])

class TestTemporaryDirectory(unittest.TestCase):
    def test_temporary_directory(self):
        with TemporaryDirectory() as tmp_dir:
            self.assertTrue(os.path.isdir(tmp_dir))

        self.assertFalse(os.path.isdir(tmp_dir))
        self.assertFalse(os.path.exists(tmp_dir))

class TestCache(TestWithData):
    def test_basics(self):
        cache = Cache()

        self.assertEqual(len(cache), 0)

        data = np.random.randn(10, 10, 3)
        cache['foo'] = data

        self.assertEqual(list(cache.keys()), ['foo'])
        np.testing.assert_allclose(cache['foo'], data)

        del cache['foo']

        self.assertEqual(len(cache), 0)

    def test_cache_to_disk(self):
        cache = Cache(cache_location = self.data_path)

        self.assertEqual(len(cache), 0)

        data = np.random.randn(10, 10, 3)
        cache['foo'] = data

        self.assertEqual(glob.glob(os.path.join(self.data_path, '*.*')),
                         [os.path.join(self.data_path, 'foo.npy')])

        self.assertEqual(list(cache.keys()), ['foo'])
        np.testing.assert_allclose(cache['foo'], data)

        cache = Cache(cache_location = self.data_path)
        self.assertEqual(cache._cache, {})
        self.assertEqual(list(cache.keys()), ['foo'])
        np.testing.assert_allclose(cache['foo'], data)

        del cache['foo']
        self.assertEqual(len(cache), 0)
        self.assertEqual(glob.glob(os.path.join(self.data_path, '*.*')),
                         [])

    def test_cache_to_disk_nonexisting_location(self):
        cache_location = os.path.join(self.data_path, 'cache')
        cache = Cache(cache_location = cache_location)

        self.assertEqual(len(cache), 0)

        data = np.random.randn(10, 10, 3)
        cache['foo'] = data

        self.assertEqual(glob.glob(os.path.join(cache_location, '*.*')),
                         [os.path.join(cache_location, 'foo.npy')])

        self.assertEqual(list(cache.keys()), ['foo'])
        np.testing.assert_allclose(cache['foo'], data)

        cache = Cache(cache_location = cache_location)
        self.assertEqual(cache._cache, {})
        self.assertEqual(list(cache.keys()), ['foo'])
        np.testing.assert_allclose(cache['foo'], data)

        del cache['foo']
        self.assertEqual(len(cache), 0)
        self.assertEqual(glob.glob(os.path.join(cache_location, '*.*')),
                         [])

    def test_pickle_cache(self):
        cache = Cache()

        self.assertEqual(len(cache), 0)

        data = np.random.randn(10, 10, 3)
        cache['foo'] = data

        self.assertEqual(list(cache.keys()), ['foo'])
        np.testing.assert_allclose(cache['foo'], data)

        cache2 = self.pickle_and_reload(cache)
        self.assertEqual(cache2._cache, {})
        self.assertEqual(len(cache2), 0)

    def test_pickle_cache_with_location(self):
        cache = Cache(cache_location = self.data_path)

        self.assertEqual(len(cache), 0)

        data = np.random.randn(10, 10, 3)
        cache['foo'] = data

        self.assertEqual(glob.glob(os.path.join(self.data_path, '*.*')),
                         [os.path.join(self.data_path, 'foo.npy')])

        self.assertEqual(list(cache.keys()), ['foo'])
        np.testing.assert_allclose(cache['foo'], data)

        cache2 = self.pickle_and_reload(cache)
        self.assertEqual(cache2._cache, {})
        self.assertEqual(len(cache2), 1)
        np.testing.assert_allclose(cache2['foo'], data)

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

ext: py

max_stars_repo_name: BigHeartDB/pysaliency

max_issues_repo_name: BigHeartDB/pysaliency

max_forks_repo_name: BigHeartDB/pysaliency

reason: import numpy


text:
# train for lines

import argparse
import numpy as np
import torch
import torch.nn as nn
import torch.nn.functional as F
import torch.optim as optim
from torchvision import datasets, transforms as tf
from torch.utils.data import DataLoader
from torch.utils.data import Dataset
import os, glob
import json
import cv2
import torchvision
from torchvision.io import read_image, ImageReadMode
from torch.optim.lr_scheduler import StepLR
import pathlib
from torchsummary import summary
from torchvision.models import resnet
import time
from scipy import ndimage,io
from skimage import io, transform

from model_lines import dla34
from torchvision import transforms
from torch.utils.data import DataLoader
from datasets import ChartDataset,collate
from loss import FocalLoss, RegL1Loss

root_dir_synth="/content/gdrive/MyDrive/CV Final/ICPR_ChartCompetition2020_AdobeData"
root_dir_pmc = "/content/gdrive/MyDrive/CV Final/ICPR2020_CHARTINFO_UB_PMC_TRAIN_v1.21/"
chartype="line"
def train(root_dir=root_dir_synth,dataset="synth", chart_type=None, img_size=(1024, 1024),path=None,epochs=5):
    model_path=path
    trainset = ChartDataset(root_dir=None, dataset=None, chart_type=chart_type, img_size=(1024, 1024),
                            heatmap_size=(256, 256))
    trainloader = torch.utils.data.DataLoader(trainset, batch_size=3, shuffle=False, collate_fn=collate, num_workers=2)
    model=dla34()
    criterion_hm = FocalLoss()
    criterion_wh = RegL1Loss()
    #criterion_wh = nn.MSELoss()
    criterion_reg = RegL1Loss()
    optimizer = optim.Adam(model.parameters(),lr=0.0001)
    loss_meter, it = 0, 0
    for i in range(epochs):
        print("epoch:: ", i)
    for j, data in enumerate(trainloader):
        model.train()
        # if torch.cuda.is_available():
        # img,heatmaps,points,boxes,wh,reg
        input_img = data[0].to(device)
        input_hm = data[1].to(device)
        input_wh = data[3].to(device)
        input_reg = data[4].to(device)
        input_regmask = data[5].to(device)
        input_ind = data[6].to(device)
        model.to(device)

        optimizer.zero_grad()

        # print(target.shape)
        output_hm, output_hw, output_reg = model(input_img)
        # print(output_hm, output_hw, output_reg)
        hm_loss = criterion_hm(output_hm, input_hm)
        # wh_loss = criterion_wh(input_wh, output_hw, input_regmask, input_ind)
        # reg_loss = criterion_reg(input_reg, output_reg,input_regmask, input_ind)
        # loss = my_loss(output, target)

        print("Loss: ", hm_loss)

        hm_loss.backward()
        optimizer.step()
        torch.cuda.empty_cache()
        loss_meter += hm_loss.item()

    model_save_state = {}
    path = model_path
    model_save_state = {
        'model_state_dict': model.state_dict(),
        'optimizer_state_dict': optimizer.state_dict(),
        'loss': hm_loss}
    torch.save(model_save_state, path)

reason: import numpy,from scipy


text:
# -*- coding: utf-8 -*-
"""
Created on Thu Oct  1 11:22:35 2015

@author: jmmauricio
"""

import numpy as np
import xlrd
import pandas as pd

def losses(i_rms, m, fp, T_a, params):

    a_i = params['a_i']
    b_i = params['b_i']
    c_i = params['c_i']
    d_i = params['d_i']
    e_i = params['e_i']
    a_d = params['a_d']
    b_d = params['b_d']
    c_d = params['c_d']
    d_d = params['d_d']
    e_d = params['e_d']

    R_th_igbt =params['R_th_igbt']
    R_th_diode=params['R_th_diode']
    R_th_igbt_case=params['R_th_igbt_case']
    R_th_diode_case=params['R_th_diode_case']
    R_th_sink =params['R_th_sink']

    p_igbt  = a_i + (b_i - c_i*m*fp)*i_rms + (d_i - e_i*m*fp)*i_rms*i_rms
    p_diode = a_d + (b_d - c_d*m*fp)*i_rms + (d_d - e_d*m*fp)*i_rms*i_rms

    p_switch = p_igbt + p_diode

    T_sink       = T_a          + p_switch*R_th_sink;
    T_case_igbt  = T_sink       + p_igbt *R_th_igbt_case;
    T_case_diode = T_sink       + p_diode*R_th_diode_case;
    T_igbt       = T_case_igbt  + p_igbt *R_th_igbt;
    T_diode      = T_case_diode + p_diode*R_th_diode;

    #print(R_th_igbt,R_th_igbt_case,T_case_igbt,T_sink)

    powers = dict(p_igbt=p_igbt,
                  p_diode=p_diode)

    temperatures = dict(T_igbt=T_igbt,
                        T_diode=T_diode,
                        T_case_igbt=T_case_igbt,
                        T_case_diode=T_case_diode,
                        T_sink=T_sink,
                        T_igbt_deg=T_igbt-273.15,
                        T_diode_deg=T_diode-273.15,
                        T_sink_deg=T_sink-273.15)

    return powers, temperatures

def vscthmodel(i_rms, m, fp, T_a, params):

    a_i = params['a_i']
    b_i = params['b_i']
    c_i = params['c_i']
    d_i = params['d_i']
    e_i = params['e_i']
    a_d = params['a_d']
    b_d = params['b_d']
    c_d = params['c_d']
    d_d = params['d_d']
    e_d = params['e_d']

    R_th_igbt_sink=params['R_th_igbt_sink']
    R_th_diode_sink=params['R_th_diode_sink']
    R_th_sink_a =params['R_th_sink_a']

    p_igbt  = a_i + (b_i - c_i*m*fp)*i_rms + (d_i - e_i*m*fp)*i_rms*i_rms
    p_diode = a_d + (b_d - c_d*m*fp)*i_rms + (d_d - e_d*m*fp)*i_rms*i_rms

    p_switch = p_igbt + p_diode

    T_sink       = T_a          + p_switch*R_th_sink_a
    T_igbt  = T_sink       + p_igbt *R_th_igbt_sink
    T_diode = T_sink       + p_diode*R_th_diode_sink

    #print(R_th_igbt,R_th_igbt_case,T_case_igbt,T_sink)

    powers = dict(p_igbt=p_igbt,
                  p_diode=p_diode)

    temperatures = dict(T_igbt=T_igbt,
                        T_diode=T_diode,
                        T_sink=T_sink,
                        T_igbt_deg=T_igbt ,
                        T_diode_deg=T_diode ,
                        T_sink_deg=T_sink)

    return powers, temperatures

def man2param_1(man_electric, man_thermal):
    '''

    Input
    -----

    List of 5 tuples

   [
   [i_1,m_1,cosphi_1,p_igbt_1,p_diode_1],
   [i_2,m_1,cosphi_1,p_igbt_2,p_diode_2],
   [i_3,m_1,cosphi_1,p_igbt_3,p_diode_3],
   [i_4,m_2,cosphi_2,p_igbt_4,p_diode_4],
   [i_5,m_2,cosphi_2,p_igbt_5,p_diode_5]
   ]

    '''

    k2deg = 273.15

    i_1 = man_electric[0][0]
    i_2 = man_electric[1][0]
    i_3 = man_electric[2][0]
    i_4 = man_electric[3][0]
    i_5 = man_electric[4][0]

    p_igbt_1 = man_electric[0][3]
    p_igbt_2 = man_electric[1][3]
    p_igbt_3 = man_electric[2][3]
    p_igbt_4 = man_electric[3][3]
    p_igbt_5 = man_electric[4][3]

    p_diode_1 = man_electric[0][4]
    p_diode_2 = man_electric[1][4]
    p_diode_3 = man_electric[2][4]
    p_diode_4 = man_electric[3][4]
    p_diode_5 = man_electric[4][4]

    m_1 = man_electric[0][1]
    m_2 = man_electric[1][1]
    m_3 = man_electric[2][1]
    m_4 = man_electric[3][1]
    m_5 = man_electric[4][1]

    cosphi_1 = man_electric[0][2]
    cosphi_2 = man_electric[1][2]
    cosphi_3 = man_electric[2][2]
    cosphi_4 = man_electric[3][2]
    cosphi_5 = man_electric[4][2]

    alpha_1 = m_1*cosphi_1
    alpha_2 = m_2*cosphi_2
    alpha_3 = m_3*cosphi_3
    alpha_4 = m_4*cosphi_4
    alpha_5 = m_5*cosphi_5

    A = np.array(
    [
    [1,  i_1,  -i_1*alpha_1,  i_1**2, -i_1**2*alpha_1],
    [1,  i_2,  -i_2*alpha_2,  i_2**2, -i_2**2*alpha_2],
    [1,  i_3,  -i_3*alpha_3,  i_3**2, -i_3**2*alpha_3],
    [1,  i_4,  -i_4*alpha_4,  i_4**2, -i_4**2*alpha_4],
    [1,  i_5,  -i_5*alpha_5,  i_5**2, -i_5**2*alpha_5]
    ]
    )

    print(A)

    b_igbt = np.array(
    [
    [p_igbt_1],
    [p_igbt_2],
    [p_igbt_3],
    [p_igbt_4],
    [p_igbt_5]
    ]
    )

    x = np.linalg.solve(A, b_igbt)

    a_i  = x[0]
    b_i = x[1]
    c_i = x[2]
    d_i = x[3]
    e_i = x[4]

    b_diode = np.array(
    [
    [p_diode_1],
    [p_diode_2],
    [p_diode_3],
    [p_diode_4],
    [p_diode_5]
    ]
    )

    x = np.linalg.solve(A, b_diode)

    a_d   = x[0]
    b_d = x[1]
    c_d = x[2]
    d_d = x[3]
    e_d = x[4]

    points=[
            [1, i_1,m_1,cosphi_1,alpha_1,p_igbt_1,p_diode_1],
            [2, i_2,m_2,cosphi_2,alpha_2,p_igbt_2,p_diode_2],
            [3, i_3,m_3,cosphi_3,alpha_3,p_igbt_3,p_diode_3],
            [4, i_4,m_4,cosphi_4,alpha_4,p_igbt_4,p_diode_4],
            [5, i_5,m_5,cosphi_5,alpha_5,p_igbt_5,p_diode_5]
            ]

    # man_thermal += [[p_igbt,p_diode,T_igbt+k2deg,T_diode+k2deg,T_sink+k2deg,T_a+k2deg]]
    idx = 0
    p_igbt   = man_thermal[idx][0]
    p_diode  = man_thermal[idx][1]

    T_igbt   = man_thermal[idx][2]
    T_diode  = man_thermal[idx][3]
    T_sink   = man_thermal[idx][4]
    T_a      = man_thermal[idx][5]

    p_switch = p_igbt + p_diode

    R_th_igbt_sink = (T_igbt-T_sink)/p_igbt

    R_th_sink_a = (T_sink-(T_a))/p_switch

    idx = 0
    p_diode = man_thermal[idx][1]
    T_diode  = man_thermal[idx][3]
    T_sink   = man_thermal[idx][4]
    R_th_diode_sink = (T_diode-T_sink)/p_diode

#    print(tabulate(points,tablefmt='latex'))

    params = dict(
    a_i = a_i,
    b_i = b_i,
    c_i = c_i,
    d_i = d_i,
    e_i = e_i,
    a_d = a_d,
    b_d = b_d,
    c_d = c_d,
    d_d = d_d,
    e_d = e_d,
    R_th_igbt_sink = R_th_igbt_sink,
    R_th_diode_sink = R_th_diode_sink,
    R_th_sink_a = R_th_sink_a,
    )

    return params

def man2param(man_data, validation=False, method='lsq1'):
    '''

    Input
    -----

    List of dictionaries, each dictionary with the following information:

     {'i': 29.0,'m':1.0,'cosphi':1.0,'p_i': 27.0,'p_d':6.99,'T_i': 54.0,'T_d': 51.0, 'T_c': 49.0, 'T_s':46.0, 'T_a':40.0 }

    'i' : RMS current (A)
    'm' : Modulator peak value (pu)
    'cosphi' : Power factor (-)
    'p_i' : IGBT power loss (W)
    'p_d' : Diode power loss (W)
    'T_i' : IGBT junture temperature (Celcius degrees)
    'T_d' : Diode junture temperature (Celcius degrees)
    'T_c' : Case temperature (Celcius degrees)
    'T_s' : Heatsink temperature (Celcius degrees)
    'T_a' : Ambient temperature (Celcius degrees)
    'Tau_h' : Sink temperature time constant (s)
    '''

    k2deg = 273.15

    N = len(man_data)
    A = np.zeros((N,5))
    b_i = np.zeros((N,1))
    b_d = np.zeros((N,1))

    A_th = np.zeros((N,5))

    R_th_igbt_sink_k = 0.0
    R_th_sink_a_k = 0.0
    R_th_diode_sink_k = 0.0

    k = 0
    k_th = 0.0
    for item in man_data:

        i_k = item['i']
        m_k = item['m']
        cosphi_k = item['cosphi']
        alpha_k = m_k*cosphi_k

        p_i_k = item['p_i']
        p_d_k = item['p_d']

        A[k,:] = np.array([1,  i_k,  -i_k*alpha_k,  i_k**2, -i_k**2*alpha_k])
        b_i[k,:] = np.array([p_i_k])
        b_d[k,:] = np.array([p_d_k])

        k += 1

        if 'T_i' in item:
            T_i = item['T_i']
            T_d = item['T_d']
            T_c = item['T_c']
            T_s = item['T_s']
            T_a = item['T_a']

            p_switch_k = p_i_k + p_d_k

            R_th_igbt_sink_k += (T_i-T_s)/p_i_k
            R_th_sink_a_k += (T_s-(T_a))/p_switch_k
            R_th_diode_sink_k += (T_d-T_s)/p_d_k

            k_th += 1.0

    R_th_igbt_sink = R_th_igbt_sink_k/k_th
    R_th_sink_a = R_th_sink_a_k/k_th
    R_th_diode_sink = R_th_diode_sink_k/k_th

   # x_i = np.linalg.solve(A, b_i)
   # x_d = np.linalg.solve(A, b_d)

    if method == 'lsq1':
        x_i = np.linalg.lstsq(A, b_i, rcond=None)[0]
        x_d = np.linalg.lstsq(A, b_d, rcond=None)[0]

    if method == 'lsq2':
        b = np.hstack((b_i,b_d))
        x = np.linalg.lstsq(A, b, rcond=None)[0]
        print(b)
        print(x)
        x_i = x[0:5,0]
        x_d = x[0:5,1]

    params = dict(
    a_i = x_i[0],
    b_i = x_i[1],
    c_i = x_i[2],
    d_i = x_i[3],
    e_i = x_i[4],
    a_d = x_d[0],
    b_d = x_d[1],
    c_d = x_d[2],
    d_d = x_d[3],
    e_d = x_d[4],
    R_th_igbt_sink = R_th_igbt_sink,
    R_th_diode_sink = R_th_diode_sink,
    R_th_sink_a = R_th_sink_a,
    )

    if validation:
        for item in man_data:
            i_rms = item['i']
            m = item['m']
            cosphi = item['cosphi']
            T_a = item['T_a']
            powers, temperatures  = vscthmodel(i_rms, m, cosphi, T_a, params)
            p_i_k, p_d_k, T_i_k, T_d_k, T_s_k =  float(item['p_i']),float(item['p_d']),float(item['T_i']),float(item['T_d']),float(item['T_s'])
            p_i, p_d, T_i, T_d, T_s =  float(powers['p_igbt']),float(powers['p_diode']),float(temperatures['T_igbt']),float(temperatures['T_diode']),float(temperatures['T_sink'])
            print(f'{p_i_k:0.2f}, {p_d_k:0.2f}, {T_i_k:0.2f}, {T_d_k:0.2f}, {T_s_k:0.2f}')
            print(f'{p_i:0.2f}, {p_d:0.2f}, {T_i:0.2f}, {T_d:0.2f}, {T_s:0.2f}')
            eps_pi = (p_i_k - p_i)/p_i_k*100.0
            eps_pd = (p_d_k - p_d)/p_d_k*100.0
            eps_ti = (T_i_k - T_i)/T_i_k*100.0
            eps_td = (T_d_k - T_d)/T_d_k*100.0
            eps_ts = (T_s_k - T_s)/T_s_k*100.0
            print(f'{eps_pi:<0.2f}%, {eps_pd:0.2f}%, {eps_ti:0.2f}%, {eps_td:0.2f}%, {eps_ts:0.2f}%')
            print(f'-------------------------------------------------------------------------------')
    return params

    # powers = dict(p_igbt=p_igbt,
    #               p_diode=p_diode)

    # temperatures = dict(T_igbt=T_igbt,
    #                     T_diode=T_diode,
    #                     T_sink=T_sink,
    #                     T_igbt_deg=T_igbt ,
    #                     T_diode_deg=T_diode ,
    #                     T_sink_deg=T_sink)

def semisel_xls(file,shs_for_param, shs_for_validate, T_a_deg):
    k2deg = 273.15

    wb = xlrd.open_workbook(file)
    man_electric = []
    man_thermal  = []

    for sh_num in shs_for_param:
        sh = wb.sheet_by_index(sh_num)

        V_dc  = float(sh.cell_value(0,1).split(' ')[0])
        V_ac  = float(sh.cell_value(1,1).split(' ')[0])
        i_rms = float(sh.cell_value(2,1).split(' ')[0])
        freq  = float(sh.cell_value(4,1).split(' ')[0])
        fp    = float(sh.cell_value(5,1))
        f_sw  = float(sh.cell_value(7,1).split(' ')[0])

        p_igbt  = float(sh.cell_value(15,1).split(' ')[0])
        p_diode = float(sh.cell_value(18,1).split(' ')[0])
        T_igbt  = float(sh.cell_value(23,1).split(' ')[0])
        T_diode = float(sh.cell_value(24,1).split(' ')[0])
        T_sink= float(sh.cell_value(21,1).split(' ')[0])

        m = V_ac*np.sqrt(2)/V_dc

        man_electric += [[i_rms,m,fp,p_igbt,p_diode]]
        man_thermal += [[p_igbt,p_diode,T_igbt+k2deg,T_diode+k2deg,T_sink+k2deg,T_a_deg+k2deg]]

        print('{} & {} & {} & {} & {} & {} & {}'.format(i_rms, fp, p_igbt, p_diode, T_igbt, T_diode, T_sink))

    params = man2param(man_electric,man_thermal)

    print('\midrule')
    for sh_num in shs_for_validate:
        sh = wb.sheet_by_index(sh_num)

        V_dc  = float(sh.cell_value(0,1).split(' ')[0])
        V_ac  = float(sh.cell_value(1,1).split(' ')[0])
        i_rms = float(sh.cell_value(2,1).split(' ')[0])
        freq  = float(sh.cell_value(4,1).split(' ')[0])
        fp    = float(sh.cell_value(5,1))
        f_sw  = float(sh.cell_value(7,1).split(' ')[0])

        p_igbt  = float(sh.cell_value(15,1).split(' ')[0])
        p_diode = float(sh.cell_value(18,1).split(' ')[0])
        T_igbt  = float(sh.cell_value(23,1).split(' ')[0])
        T_diode = float(sh.cell_value(24,1).split(' ')[0])
        T_sink= float(sh.cell_value(21,1).split(' ')[0])

        m = V_ac*np.sqrt(2)/V_dc

        pows, temps    = vscthmodel(i_rms, m, fp, T_a_deg, params)

        print('{:2.1f} & {:2.2f} & {:2.1f} & {:2.1f}  & {:2.1f} & {:2.1f}   & {:2.1f} & {:2.1f}  & {:2.1f} & {:2.1f}   & {:2.1f} & {:2.1f} \\\\  '.format(i_rms,
                                                                                fp,
                                                                                p_igbt, pows['p_igbt'][0],
                                                                                p_diode, pows['p_diode'][0],
                                                                                T_igbt, temps['T_igbt_deg'][0],
                                                                                T_diode, temps['T_diode_deg'][0],
                                                                                T_sink, temps['T_sink_deg'][0],

))

    return params

def imposim_xls(file):
    k2deg = 273.15

    wb = xlrd.open_workbook(file)
    man_electric = []
    man_thermal  = []

    sh = wb.sheet_by_index(0)

    V_dc  = float(sh.cell_value(20,3))
    m     = float(sh.cell_value(26,3))

    idx = 9

    T_a_deg = float(sh.cell_value(idx+53,9))

    i_rms   = float(sh.cell_value(idx,0))
    m       = float(sh.cell_value(26,3))
    fp      = float(sh.cell_value(27,3))
    p_igbt  = float(sh.cell_value(idx,6))+ 0.5*float(sh.cell_value(idx,13))
    p_diode = float(sh.cell_value(idx,11)) + 0.5*float(sh.cell_value(idx,13))
    T_igbt  = float(sh.cell_value(idx+34,10))
    T_diode = float(sh.cell_value(idx+52,10))
    T_sink = float(sh.cell_value(idx+52,8)) + float(sh.cell_value(idx+53,9))
    man_electric += [[i_rms,m,fp,p_igbt,p_diode]]
    man_thermal += [[p_igbt,p_diode,T_igbt+k2deg,T_diode+k2deg,T_sink+k2deg,T_a_deg+k2deg]]
    print('{:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} \\\\'.format(i_rms, fp, p_igbt, p_diode, T_igbt, T_diode, T_sink))

    idx = 13
    i_rms   = float(sh.cell_value(idx,0))
    m       = float(sh.cell_value(26,3))
    fp      = float(sh.cell_value(27,3))
    p_igbt  = float(sh.cell_value(idx,6))+ 0.5*float(sh.cell_value(idx,13))
    p_diode = float(sh.cell_value(idx,11))+ 0.5*float(sh.cell_value(idx,13))
    T_igbt  = float(sh.cell_value(idx+34,10))
    T_diode = float(sh.cell_value(idx+52,10))
    T_sink = float(sh.cell_value(idx+52,8)) + float(sh.cell_value(idx+53,9))
    man_electric += [[i_rms,m,fp,p_igbt,p_diode]]
    man_thermal += [[p_igbt,p_diode,T_igbt+k2deg,T_diode+k2deg,T_sink+k2deg,T_a_deg+k2deg]]
    print('{:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} \\\\'.format(i_rms, fp, p_igbt, p_diode, T_igbt, T_diode, T_sink))

    idx = 16
    i_rms   = float(sh.cell_value(idx,0))
    m       = float(sh.cell_value(26,3))
    fp      = float(sh.cell_value(27,3))
    p_igbt  = float(sh.cell_value(idx,6)) + 0.5*float(sh.cell_value(idx,13))
    p_diode = float(sh.cell_value(idx,11))  + 0.5*float(sh.cell_value(idx,13))
    T_igbt  = float(sh.cell_value(idx+34,10))
    T_diode = float(sh.cell_value(idx+52,10))
    T_sink = float(sh.cell_value(idx+52,8)) + float(sh.cell_value(idx+53,9))
    man_electric += [[i_rms,m,fp,p_igbt,p_diode]]
    man_thermal += [[p_igbt,p_diode,T_igbt+k2deg,T_diode+k2deg,T_sink+k2deg,T_a_deg+k2deg]]
    print('{:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} \\\\'.format(i_rms, fp, p_igbt, p_diode, T_igbt, T_diode, T_sink))

    sh = wb.sheet_by_index(1)
    m     = float(sh.cell_value(26,3))

    idx = 13
    i_rms   = float(sh.cell_value(idx,0))
    m       = float(sh.cell_value(26,3))
    fp      = float(sh.cell_value(27,3))
    p_igbt  = float(sh.cell_value(idx,6))    + 0.5*float(sh.cell_value(idx,13))
    p_diode = float(sh.cell_value(idx,11))   + 0.5*float(sh.cell_value(idx,13))
    T_igbt  = float(sh.cell_value(idx+34,10))
    T_diode = float(sh.cell_value(idx+52,10))
    T_sink = float(sh.cell_value(idx+52,8)) + float(sh.cell_value(idx+53,9))
    man_electric += [[i_rms,m,fp,p_igbt,p_diode]]
    man_thermal += [[p_igbt,p_diode,T_igbt+k2deg,T_diode+k2deg,T_sink+k2deg,T_a_deg+k2deg]]
    print('{:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} \\\\'.format(i_rms, fp, p_igbt, p_diode, T_igbt, T_diode, T_sink))

    idx = 16
    i_rms   = float(sh.cell_value(idx,0))
    m       = float(sh.cell_value(26,3))
    fp      = float(sh.cell_value(27,3))
    p_igbt  = float(sh.cell_value(idx,6)) + 0.5*float(sh.cell_value(idx,13))
    p_diode = float(sh.cell_value(idx,11)) + 0.5*float(sh.cell_value(idx,13))
    T_igbt  = float(sh.cell_value(idx+34,10))
    T_diode = float(sh.cell_value(idx+52,10))
    T_sink = float(sh.cell_value(idx+52,8)) + float(sh.cell_value(idx+53,9))
    man_electric += [[i_rms,m,fp,p_igbt,p_diode]]
    man_thermal += [[p_igbt,p_diode,T_igbt+k2deg,T_diode+k2deg,T_sink+k2deg,T_a_deg+k2deg]]
    print('{:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} & {:2.1f} \\\\'.format(i_rms, fp, p_igbt, p_diode, T_igbt, T_diode, T_sink))

    print(np.array(man_electric))
    params = man2param(man_electric,man_thermal)

def imposim2xls(file):
    k2deg = 273.15

    df_1 = pd.read_excel(file, sheetname=0, skiprows=8)
    df_2 = pd.read_excel(file, sheetname=1, skiprows=8)
    print(df_1)

    idx1 = 3
    idx2 = 10
    idx3 = 16

    man_electric = [
                   [df_1.i_rms[idx1],df_1.m[idx1],df_1.fp[idx1],df_1.p_igbt[idx1],df_1.p_diode[idx1]],
                   [df_1.i_rms[idx2],df_1.m[idx2],df_1.fp[idx2],df_1.p_igbt[idx2],df_1.p_diode[idx2]],
                   [df_1.i_rms[idx3],df_1.m[idx3],df_1.fp[idx3],df_1.p_igbt[idx3],df_1.p_diode[idx3]],
                   [df_2.i_rms[idx2],df_2.m[idx2],df_2.fp[idx2],df_2.p_igbt[idx2],df_2.p_diode[idx2]],
                   [df_2.i_rms[idx3],df_2.m[idx3],df_2.fp[idx3],df_2.p_igbt[idx3],df_2.p_diode[idx3]],
                   ]

    man_thermal = [
                  [df_1.p_igbt[idx2],df_1.p_diode[idx2],df_1.T_igbt[idx2]+k2deg,df_1.T_diode[idx2]+k2deg,df_1.T_sink[idx2]+k2deg,df_1.T_a[idx2]+k2deg],
                  ]

    params =  man2param(man_electric, man_thermal)

    return params

if __name__ == '__main__':
    k2deg = 273.15
    file = '/home/jmmauricio/Documents/public/jmmauricio6/RESEARCH/ARTICLES/doing/vsc_model/code/semikron_100kva/semikron_SKiiP38GB12E4V1.xls'
#    shs_for_param = [0,2,4,5,6]
#    shs_for_validate = [7,8,9,10]
#    T_a_deg = 40.0+k2deg
#    params = semisel_xls(file,shs_for_param, shs_for_validate, T_a_deg)
#
#    print(params)
    file = '/home/jmmauricio/Documents/public/jmmauricio6/RESEARCH/ARTICLES/doing/vsc_model/code/imposim/imposim_F800R17.xls'
    params = imposim2xls(file)
    print(params)

    T_a = 40+ k2deg
    R_th_diode_sink = params['R_th_diode_sink']
    R_th_igbt_sink = params['R_th_igbt_sink']
    R_th_sink_a = params['R_th_sink_a']
    a_d = params['a_d']
    a_i = params['a_i']
    b_d = params['b_d']
    b_i = params['b_i']
    c_d = params['c_d']
    c_i = params['c_i']
    d_d = params['d_d']
    d_i = params['d_i']
    e_d = params['e_d']
    e_i = params['e_i']

    i_rms = 520

    m = 400*np.sqrt(2)/800
    fp =  1.0
    p_igbt  = 1.000*(a_i + (b_i - c_i*m*fp)*i_rms + (d_i - e_i*m*fp)*i_rms*i_rms)
    p_diode = 1.000*(a_d + (b_d - c_d*m*fp)*i_rms + (d_d - e_d*m*fp)*i_rms*i_rms)

    p_switch = p_igbt + p_diode + i_rms**2*0.000

    print('p_igbt',p_igbt)
    print('p_diode',p_diode)
    print('p_switch', p_switch)
    T_sink_0       = T_a    + p_switch*R_th_sink_a;
    T_igbt_0       = T_sink_0 + p_igbt *R_th_igbt_sink;
    T_diode_0      = T_sink_0 + p_diode*R_th_diode_sink;

    print('T_sink_0', T_sink_0-k2deg)
    print('T_igbt_0', T_igbt_0-k2deg)
    print('T_diode_0', T_diode_0-k2deg)
#
##    C_th_diode= 10
##    C_th_diode_case= 2
##    C_th_igbt= 18
##    C_th_igbt_case= 5
##    C_th_sink= 6000.0
##    R_th_diode= 0.01979045401629802
##    R_th_diode_case= 0.018
##    R_th_igbt= 0.009765625
##    R_th_igbt_case= 0.009
##    R_th_sink= 0.007
##    a_d = 143.48507451
##    a_i = 421.02132341
##    b_d = 0.589627
##    b_i = 0.55708434
##    c_d = 0.18337165
##    c_i =-0.12254324
##    d_d = 0.00026235
##    d_i = 0.00089385
##    e_d = 0.00021407
##    e_i =-0.00041411
##    T_a = 35.0+273
##    params ={'C_th_diode': C_th_diode,
##     'C_th_diode_case': C_th_diode_case,
##     'C_th_igbt': C_th_igbt,
##     'C_th_igbt_case': C_th_igbt_case,
##     'C_th_sink': C_th_sink,
##     'R_th_diode': R_th_diode,
##     'R_th_diode_case': R_th_diode_case,
##     'R_th_igbt': R_th_igbt,
##     'R_th_igbt_case':R_th_igbt_case,
##     'R_th_sink': R_th_sink,
##     'a_d':  a_d,
##     'a_i':  a_i,
##     'b_d':  b_d,
##     'b_i':  b_i,
##     'c_d':  c_d,
##     'c_i':  c_i,
##     'd_d':  d_d,
##     'd_i':  d_i,
##     'e_d':  e_d,
##     'e_i':  e_i,
##     'm':0.85,
##     'fp':0.8,
##     'i_rms':1400,
##      'T_a':35+273.3,
##            }
##
##
##    i_rms = 1500
##    fp = 0.85
##    m = 0.8
##    pows, temps = losses(i_rms, m, fp, T_a, params)
##    print(pows)
##    print(temps)

max_stars_repo_name: pydgrid/pydgrid

max_issues_repo_name: pydgrid/pydgrid

max_forks_repo_name: pydgrid/pydgrid


text:
"""Utility functions."""

from typing import Tuple

import numpy

from rio_color.operations import parse_operations
from rio_color.utils import scale_dtype, to_math_type

from rio_tiler.utils import linear_rescale, _chunks

def _postprocess(
    tile: numpy.ndarray,
    mask: numpy.ndarray,
    rescale: str = None,
    color_formula: str = None,
) -> Tuple[numpy.ndarray, numpy.ndarray]:
    if rescale:
        rescale_arr = list(map(float, rescale.split(",")))
        rescale_arr = list(_chunks(rescale_arr, 2))
        if len(rescale_arr) != tile.shape[0]:
            rescale_arr = ((rescale_arr[0]),) * tile.shape[0]
        for bdx in range(tile.shape[0]):
            tile[bdx] = numpy.where(
                mask,
                linear_rescale(
                    tile[bdx], in_range=rescale_arr[bdx], out_range=[0, 255]
                ),
                0,
            )
        tile = tile.astype(numpy.uint8)

    if color_formula:
        # make sure one last time we don't have
        # negative value before applying color formula
        tile[tile < 0] = 0
        for ops in parse_operations(color_formula):
            tile = scale_dtype(ops(to_math_type(tile)), numpy.uint8)

    return tile, mask

def intensity_range(image, range_values='image', clip_negative=False):
    """Return image intensity range (min, max) based on desired value type.
    Parameters
    ----------
    image : array
        Input image.
    range_values : str or 2-tuple, optional
        The image intensity range is configured by this parameter.
        The possible values for this parameter are enumerated below.
        'image'
            Return image min/max as the range.
        'dtype'
            Return min/max of the image's dtype as the range.
        dtype-name
            Return intensity range based on desired `dtype`. Must be valid key
            in `DTYPE_RANGE`. Note: `image` is ignored for this range type.
        2-tuple
            Return `range_values` as min/max intensities. Note that there's no
            reason to use this function if you just want to specify the
            intensity range explicitly. This option is included for functions
            that use `intensity_range` to support all desired range types.
    clip_negative : bool, optional
        If True, clip the negative range (i.e. return 0 for min intensity)
        even if the image dtype allows negative values.
    """
    _integer_types = (numpy.byte, numpy.ubyte,      # 8 bits
                  numpy.short, numpy.ushort,        # 16 bits
                  numpy.intc, numpy.uintc,          # 16 or 32 or 64 bits
                  numpy.int_, numpy.uint,           # 32 or 64 bits
                  numpy.longlong, numpy.ulonglong)  # 64 bits
    _integer_ranges = {t: (numpy.iinfo(t).min, numpy.iinfo(t).max)
                   for t in _integer_types}
    dtype_range = {numpy.bool_: (False, True),
               numpy.bool8: (False, True),
               numpy.float16: (-1, 1),
               numpy.float32: (-1, 1),
               numpy.float64: (-1, 1)}
    dtype_range.update(_integer_ranges)

    DTYPE_RANGE = dtype_range.copy()
    DTYPE_RANGE.update((d.__name__, limits) for d, limits in dtype_range.items())
    DTYPE_RANGE.update({'uint10': (0, 2 ** 10 - 1),
                        'uint12': (0, 2 ** 12 - 1),
                        'uint14': (0, 2 ** 14 - 1),
                        'bool': dtype_range[numpy.bool_],
                        'float': dtype_range[numpy.float64]})

    if range_values == 'dtype':
        range_values = image.dtype.type

    if range_values == 'image':
        i_min = numpy.min(image)
        i_max = numpy.max(image)
    elif range_values in DTYPE_RANGE:
        i_min, i_max = DTYPE_RANGE[range_values]
        if clip_negative:
            i_min = 0
    else:
        i_min, i_max = range_values
    return i_min, i_max

def rescale_intensity(image, in_range='image', out_range='dtype'):
    """Return image after stretching or shrinking its intensity levels.
    The desired intensity range of the input and output, `in_range` and
    `out_range` respectively, are used to stretch or shrink the intensity range
    of the input image. See examples below.
    Parameters
    ----------
    image : array
        Image array.
    in_range, out_range : str or 2-tuple, optional
        Min and max intensity values of input and output image.
        The possible values for this parameter are enumerated below.
        'image'
            Use image min/max as the intensity range.
        'dtype'
            Use min/max of the image's dtype as the intensity range.
        dtype-name
            Use intensity range based on desired `dtype`. Must be valid key
            in `DTYPE_RANGE`.
        2-tuple
            Use `range_values` as explicit min/max intensities.
    Returns
    -------
    out : array
        Image array after rescaling its intensity. This image is the same dtype
        as the input image.
    See Also
    --------
    equalize_hist
    Examples
    --------
    By default, the min/max intensities of the input image are stretched to
    the limits allowed by the image's dtype, since `in_range` defaults to
    'image' and `out_range` defaults to 'dtype':
    >>> image = np.array([51, 102, 153], dtype=np.uint8)
    >>> rescale_intensity(image)
    array([  0, 127, 255], dtype=uint8)
    It's easy to accidentally convert an image dtype from uint8 to float:
    >>> 1.0 * image
    array([ 51., 102., 153.])
    Use `rescale_intensity` to rescale to the proper range for float dtypes:
    >>> image_float = 1.0 * image
    >>> rescale_intensity(image_float)
    array([0. , 0.5, 1. ])
    To maintain the low contrast of the original, use the `in_range` parameter:
    >>> rescale_intensity(image_float, in_range=(0, 255))
    array([0.2, 0.4, 0.6])
    If the min/max value of `in_range` is more/less than the min/max image
    intensity, then the intensity levels are clipped:
    >>> rescale_intensity(image_float, in_range=(0, 102))
    array([0.5, 1. , 1. ])
    If you have an image with signed integers but want to rescale the image to
    just the positive range, use the `out_range` parameter:
    >>> image = np.array([-10, 0, 10], dtype=np.int8)
    >>> rescale_intensity(image, out_range=(0, 127))
    array([  0,  63, 127], dtype=int8)
    """
    dtype = image.dtype.type

    imin, imax = intensity_range(image, in_range)
    omin, omax = intensity_range(image, out_range, clip_negative=(imin >= 0))

    # Fast test for multiple values, operations with at least 1 NaN return NaN
    if numpy.isnan(imin + imax + omin + omax):
        return
    image = numpy.clip(image, imin, imax)

    if imin != imax:
        image = (image - imin) / float(imax - imin)
    return numpy.asarray(image * (omax - omin) + omin, dtype=dtype)

max_stars_repo_name: betonr/remotepixel-tiler

max_issues_repo_name: betonr/remotepixel-tiler

max_forks_repo_name: betonr/remotepixel-tiler


text:
import copy
import os, sys
import random
import numpy as np
import json

import tensorflow as tf
from tensorflow.keras.utils import Sequence, to_categorical
from glob import glob
from absl import flags
from tqdm import tqdm

from split import shuffle_arrays, shuffle_imu
from utils import gyro_conv
import myFlags
FLAGS = flags.FLAGS

class TotalCapture():
    Labels =['Head','Sternum','Pelvis','L_UpArm','R_UpArm',\
        'L_LowArm','R_LowArm','L_UpLeg','R_UpLeg','L_LowLeg',\
        'R_LowLeg','L_Foot','R_Foot']
    ObjectCategories = len(Labels)
    SymmeLabels = ['Head','Sternum','Pelvis','UpArm','LowArm',\
        'UpLeg','LowLeg','Foot']
    window_size = 120
    window_interval = 15
    ignore_frames = 240 # remove T-pose frames

    def __init__(self, root=None, used_imus=list(range(13)), label_imus=list(range(13))):
        self.DATASET_ROOT_PATH = FLAGS.datapath
        self.used_imus = copy.deepcopy(used_imus)
        self.label_imus = copy.deepcopy(label_imus)
        self.nb_classes = len(set(label_imus))
        self.root_joint = "Pelvis" if root is None else root

    def path(self, data_type=-1, subject_id=-1, scene=-1):
        p = self.DATASET_ROOT_PATH

        # data_type is 'orientation' or 'acceleration'
        if data_type == -1:
            return p
        p0 = os.path.join(p, data_type)

        if subject_id == -1:
            return p
        scene_name = 's'+str(subject_id)
        p = os.path.join(p0, scene_name)

        if scene == -1:
            return p
        scene_name = scene_name + '_' + scene + '.npy'
        p = os.path.join(p0, scene_name)
        return p

    def files_to_load(self, data_type, subject_ids, scenes):
        # data_type is 'orientation' or 'acceleration'
        p = os.path.join(self.DATASET_ROOT_PATH, data_type)
        dst = []
        for ins_id in subject_ids:
            for scene in scenes:
                fpath = glob(p+'/s{}_{}*.npy'.format(ins_id, scene))
                dst.extend(fpath)
        return dst

    def get_class_id(self, sub):
        return self.Labels.index(sub)

    def get_class_name(self, sub):
        if self.nb_classes == 8:
            return self.SymmeLabels[int(sub)]
        else:
            return self.Labels[int(sub)]

    def get_class_names(self, arr):
        return np.array([self.get_class_name(elem) for elem in arr])

    def get_different_modal_path(self, fpath, modal):
        oripath = os.path.dirname(fpath).split('/')[:-1] + [modal, os.path.basename(fpath)]
        return ''.join([st + '/' for st in oripath])[:-1]

class CMU_MoCap():
    Labels = ["head","thorax","lowerback","lhumerus","rhumerus",
              "lradius","rradius","lwrist","rwrist","lfemur",
              "rfemur","ltibia","rtibia","lfoot","rfoot"]
    ObjectCategories = len(Labels)
    window_size = 240
    window_interval = 30
    ignore_frames = 1

    def __init__(self, root=None, used_imus=list(range(13)), label_imus=list(range(13))):
        self.DATASET_ROOT_PATH = FLAGS.datapath
        self.used_imus = copy.deepcopy(used_imus)
        self.label_imus = copy.deepcopy(label_imus)
        self.nb_classes = len(set(label_imus))
        self.root_joint = "lowerback" if root is None else root

    def path(self, data_type=-1, subject_id=-1, scene=-1):
        p = self.DATASET_ROOT_PATH

        # data_type is 'orientation', 'acceleration', or "gyro"
        if data_type == -1:
            return p
        p0 = os.path.join(p, data_type)

        if subject_id == -1:
            return p
        scene_name = 's'+str(subject_id)
        p = os.path.join(p0, scene_name)

        if scene == -1:
            return p
        scene_name = scene_name + '_' + scene + '.npy'
        p = os.path.join(p0, scene_name)
        return p

    def files_to_load(self, data_type, subject_ids=None, scenes=None):
        # data_type is 'orientation', 'acceleration', or "gyro"
        p = os.path.join(self.DATASET_ROOT_PATH, data_type)
        dst = []
        if subject_ids in ["train_json", "valid_json", "test_json"]:
            dtype = subject_ids[:-5]
            #with open(FLAGS.scene_file, "r") as rf:
            with open(os.path.join(self.DATASET_ROOT_PATH, "walk.json"), "r") as rf:
                files = json.load(rf)
                for subj in files[dtype]:
                    dst.extend([p + "/s{:03}_{:02}.npy".format(int(subj), sce) \
                        for sce in files[dtype][subj]])
        else:
            if subject_ids is None and scenes is None:
                dst = glob(p+"/s*.npy")
            elif subject_ids is None and scenes is not None:
                for scene in scenes:
                    if type(scene) != str:
                        scene = str(scene).zfill(2)
                    fpath = glob(p+"/s*_{}.npy".format(scene))
                    dst.extend(fpath)
            elif subject_ids is not None and scenes is None:
                for ins_id in subject_ids:
                    if type(ins_id) != str:
                        ins_id = str(ins_id).zfill(3)
                    fpath = glob(p+"/s{}_*.npy".format(ins_id))
                    dst.extend(fpath)
            elif subject_ids is not None and scenes is not None:
                for ins_id in subject_ids:
                    if type(ins_id) != str:
                        ins_id = str(ins_id).zfill(3)
                    for scene in scenes:
                        if type(scene) != str:
                            scene = str(scene).zfill(2)
                        fpath = p+'/s{}_{}.npy'.format(ins_id, scene)
                        if os.path.exists(fpath):
                            dst.append(fpath)
            else:
                raise NotImplementedError("files_to_load error")
        return dst

    def get_class_id(self, sub):
        return self.Labels.index(sub)

    def get_class_name(self, sub):
        return self.Labels[int(sub)]

    def get_class_names(self, arr):
        return np.array([self.get_class_name(elem) for elem in arr])

    def get_different_modal_path(self, fpath, modal):
        oripath = os.path.dirname(fpath).split('/')[:-1] + [modal, os.path.basename(fpath)]
        return ''.join([st + '/' for st in oripath])[:-1]

class IMUsSequence(Sequence):
    def __init__(self, batch_size, dataset,
                 data_types = ['acceleration','gyro'],
                 used_imus = list(range(13)),
                 label_imus = list(range(13)),
                 random_state=None,
                 dims= 6, # sum of input dimensions
                 subject_ids=[1,2,3], # 4 for val, 5 for test
                 scenes=['freestyle','walking','acting'],
                 is_acc_norm = True):
        self.batch_size = batch_size
        self.dataset = dataset
        self.data_types = data_types
        self.used_imus = used_imus
        self.label_imus = label_imus
        self.is_acc_norm = is_acc_norm
        self.window_size = dataset.window_size
        self.window_interval = dataset.window_interval
        self.ignore_frames = dataset.ignore_frames
        self.random_state = random_state
        self.dims = dims
        self.subject_ids = copy.deepcopy(subject_ids)
        self.scenes = copy.deepcopy(scenes)

        self.nb_imus = len(self.used_imus)
        # root joint is ignored in target and label
        self.root_id = self.dataset.get_class_id(self.dataset.root_joint)
        self.nb_classes = self.get_nb_classes()
        self.nb_target_imus = len(self.used_imus) - 1

        nb_frames = []
        root_rots = []
        # load all imu data
        for dt_count, data_type in enumerate(self.data_types):
            files_to_load = self.dataset.files_to_load(data_type, self.subject_ids, self.scenes)
            files_to_load.sort()
            xs = []
            for fpath in tqdm(files_to_load, desc=data_type):
                sample = np.expand_dims(np.load(fpath), axis=-1)# (imu_id, nb_frames, quaternion, 1)
                xs.append(sample)
                if dt_count == 0:
                    nb_frames.append(sample.shape[1])
                    rot_fpath = self.dataset.get_different_modal_path(fpath, 'orientation')
                    root_rot = np.load(rot_fpath)[self.root_id] #(nb_frames, 3, 3) rot_mat
                    root_rots.append(root_rot)
            xs = np.concatenate(xs, axis=1)
            if dt_count == 0:
                stock_x = xs
            else:
                stock_x = np.concatenate([stock_x, xs], axis=2) # (nb_all_imus, frame, sum_datadims, 1)
        #root_rots = np.concatenate(xs, axis=0)
        root_rots = np.concatenate(root_rots, axis=0)
        self.nb_frames = np.asarray(nb_frames, dtype='uint16')
        if self.dataset.DATASET_ROOT_PATH.find("TotalCapture") != -1:
            self.window_first_ids =self.get_window_first_ids_removeBackWalk(files_to_load) # first frames of windows
        else:
            self.window_first_ids =self.get_window_first_ids() # first frames of windows
        self.nb_samples = len(self.window_first_ids)
        self.redefine_labels() # lift root joint to top of self.used_imus and change labels

        self.x = np.empty((self.nb_samples,) +
                self.get_x_sample_shape(), dtype=self.get_x_dtype())
        self.rot_x = np.empty((self.nb_samples,) + # rotation mat of root joint
                self.get_rot_x_sample_shape(), dtype=self.get_x_dtype())
        self.y = np.empty((self.nb_samples,) +
                self.get_y_sample_shape(), dtype=self.get_y_dtype())

        for frame_i, frame in enumerate(self.window_first_ids):
            for cat_i, cat in enumerate(self.used_imus):
                # only used IMU data is extracted
                self.x[frame_i, cat_i] = stock_x[cat, frame:frame+self.window_size] #(nb_samples, nb_used_imus, sum_datadims, 1)
                self.y[frame_i, cat_i] = self.label_imus[cat_i] #(nb_samples, nb_imus,1)
            self.rot_x[frame_i] = root_rots[frame] #(nb_samples, nb_imus,1)
        self.y = to_categorical(self.y, num_classes=self.nb_classes)

        if self.is_acc_norm: # subtract root acceleration to normalize
            self.x[:,:,:,:3] -= self.x[:,[0],:,:3]

        shuffle_arrays([self.x, self.rot_x, self.y]) # shuffle frames
        shuffle_imu([self.x, self.y], ignoreTop=True) # shuffle IMUs

    def get_x_sample_shape(self):
        return (self.nb_imus, self.window_size, self.dims, 1)

    def get_rot_x_sample_shape(self):
        return (3,3)

    def get_y_sample_shape(self):
        return (self.nb_imus,)

    def get_rot_y_sample_shape(self):
        return (self.window_size,3,3)

    def get_x_dtype(self):
        return "float32"

    def get_y_dtype(self):
        return "uint8"
        #return "float32"

    def __len__(self):
        return int(np.ceil(float(self.nb_samples) / float(self.batch_size)))

    def actual_len(self):
        # actual length without repeated augumentation
        return int(np.ceil(float(self.nb_samples) / float(self.batch_size)))

    def __getitem__(self, idx):
        idx %= self.actual_len()
        start = idx * self.batch_size

        batch_x = self.x[start:start + self.batch_size]
        ys = self.y[start:start + self.batch_size, 1:]
        rot_xs = self.rot_x[start:start + self.batch_size]
        xs = np.stack([self.augment(x, rot_x) for x, rot_x in
                zip(batch_x.astype(self.get_x_dtype()), rot_xs.astype(self.get_x_dtype()))], axis=0)
        return xs, ys

    def augment(self, x, rot_x=None):
        # x[:,:,:3] is acceleration, x[:,:,3:] is gyro
        if rot_x is not None:
            dst = np.empty_like(x)
            dst[:,:,:3] = np.matmul(rot_x.T, x[:,:,:3])
            dst[:,:,3:] = np.matmul(gyro_conv(rot_x.T), x[:,:,3:])
            return dst
        else:
            return x

    def on_epoch_end(self):
        if self.random_state is None:
            shuffle_arrays([self.x, self.rot_x, self.y]) # shuffle frames
            shuffle_imu([self.x, self.y], ignoreTop=True) # shuffle IMUs
            pass

    def get_window_first_ids(self):
        prev_frames_sum = 0
        dst = []
        for nb_frame in self.nb_frames:
            if nb_frame-self.ignore_frames*2 < self.window_size:
                #raise ValueError('number of frames is small')
                pass
            for idx in range(self.ignore_frames, nb_frame-self.ignore_frames-self.window_size, self.window_interval):
                dst.append(prev_frames_sum + idx)
            prev_frames_sum += nb_frame
        random.shuffle(dst)
        return dst

    def get_window_first_ids_removeBackWalk(self, files_to_load):
        assert len(self.nb_frames) == len(files_to_load)
        #remove_fn, remove_frame
        prev_frames_sum = 0
        tmp = []
        for nb_frame, fn in zip(self.nb_frames, files_to_load):
            if nb_frame-self.ignore_frames*2 < self.window_size:
                #raise ValueError('number of frames is small')
                pass
            for idx in range(self.ignore_frames, nb_frame-self.ignore_frames-self.window_size, self.window_interval):
                tmp.append(prev_frames_sum + idx)
            prev_frames_sum += nb_frame

        with open(os.path.join(self.dataset.DATASET_ROOT_PATH, "backwalk.json"), 'r') as f:
            remove_files = json.load(f) # {"fn": [[frame_start, frame_end],[,]..]}
        tmp_frames = 0
        remove_frames = []
        for nb_frame, fn in zip(self.nb_frames, files_to_load):
            for rm_fn in remove_files.keys():
                if fn.find(rm_fn) != -1:
                    rm_frames = remove_files[rm_fn] #[[frame_start, frame_end], [,]..]
                    for rm_frame in rm_frames:
                        rm_start = tmp_frames + rm_frame[0]
                        rm_end = tmp_frames + rm_frame[1]
                        remove_frames.append([rm_start, rm_end])
            tmp_frames += nb_frame

        dst = []
        for frame in tmp:
            should_remain = True
            for remove_frame in remove_frames:
                if frame + self.window_size > remove_frame[0] and frame < remove_frame[1]:
                    should_remain = False
            if should_remain:
                dst.append(frame)
        random.shuffle(dst)
        return dst

    def redefine_labels(self):
        root_idx = self.used_imus.index(self.root_id)
        self.used_imus.remove(self.root_id)
        self.used_imus = [self.root_id] + self.used_imus

        root_label = self.label_imus[root_idx]
        self.label_imus.remove(root_label)
        self.label_imus = np.array(self.label_imus)
        self.label_imus = [0] + list(np.where(self.label_imus > root_label,
                    self.label_imus-1, self.label_imus))

    def get_nb_classes(self):
        tmp = copy.deepcopy(self.used_imus)
        tmp.remove(self.root_id)
        return len(set(tmp))


text:
# pylint: disable=locally-disabled, multiple-statements, fixme, line-too-long, no-name-in-module, unused_import, wrong-import-order, bad-option-value
# some imports are not accepted by pylint

from models.RainbowModel import RainbowModel
import numpy as np
import matplotlib.pyplot as plt
import tensorflow as tf
import h5py
from tensorflow.keras import regularizers
from tensorflow.keras.layers import (
    Input,
    Conv2D,
    Dense,
    Flatten,
    Dropout,
    LSTM,
    GlobalMaxPooling1D,
    MaxPooling2D,
    BatchNormalization,
)
from tensorflow.keras.models import Model
from tensorflow.keras.optimizers import Adam
from sklearn.model_selection import train_test_split
from sklearn.metrics import confusion_matrix
import itertools

from datetime import datetime
import os
from utils.typing import assert_type
from utils.Window import Window
from utils.Recording import Recording
import itertools

class JensModel(RainbowModel):
    def __init__(self, **kwargs):
        """

        epochs=10
        :param kwargs:
            window_size: int
            n_features: int
            n_outputs: int
        """

        # hyper params to instance vars
        super().__init__(**kwargs)
        self.verbose = kwargs.get("verbose") or True
        self.n_epochs = kwargs.get("n_epochs") or 10
        self.learning_rate = kwargs.get("learning_rate") or 0.001
        self.model_name = "jens_model"

        self.n_features = kwargs["n_features"]
        self.n_outputs = kwargs["n_outputs"]
        self.window_size = kwargs["window_size"]

        # create model
        self.model = self._create_model()
        # Refactoring idea:
        # n_features of a neuronal net is the number of inputs, so in reality n_features = window_size * n_features
        # we could have another name for that

    def _create_model(self):

        i = Input(
            shape=(self.window_size, self.n_features, 1)
        )  # before: self.x_train[0].shape - (25, 51, 1)... before self_x_train = np.expand_dims(self.x_train[0], -1) - around the value another []
        x = Conv2D(
            32,
            (3, 3),
            strides=2,
            activation="relu",
            padding="same",
            kernel_regularizer=regularizers.l2(0.0005),
        )(i)
        x = BatchNormalization()(x)
        x = MaxPooling2D((2, 2))(x)
        x = Dropout(0.2)(x)
        x = Conv2D(
            64,
            (3, 3),
            strides=2,
            activation="relu",
            padding="same",
            kernel_regularizer=regularizers.l2(0.0005),
        )(x)
        x = BatchNormalization()(x)
        x = Dropout(0.4)(x)
        x = Conv2D(
            128,
            (3, 3),
            strides=2,
            activation="relu",
            padding="same",
            kernel_regularizer=regularizers.l2(0.0005),
        )(x)
        x = BatchNormalization()(x)
        x = MaxPooling2D((2, 2))(x)
        x = Dropout(0.2)(x)
        x = Flatten()(x)
        x = Dropout(0.2)(x)
        x = Dense(1024, activation="relu")(x)
        x = Dropout(0.2)(x)
        x = Dense(self.n_outputs, activation="softmax")(x)
        model = Model(i, x)
        model.compile(
            optimizer=Adam(learning_rate=self.learning_rate),
            loss="categorical_crossentropy",  # CategoricalCrossentropy (than we have to to the one hot encoding - to_categorical), before: "sparse_categorical_crossentropy"
            metrics=["accuracy"],
        )

        return model

    def convert(self, windows: "list[Window]") -> "tuple[np.ndarray, np.ndarray]":
        X_train, y_train = super().convert(windows)
        return np.expand_dims(X_train, -1), y_train

max_stars_repo_name: Sensors-in-Paradise/OpportunityML

max_issues_repo_name: Sensors-in-Paradise/OpportunityML

max_forks_repo_name: Sensors-in-Paradise/OpportunityML


text:
"""

Define units and operations on them.

Copyright 2012, Casey W. Stark. See LICENSE.txt for more information.

"""

from sympy import Expr, Mul, nsimplify, Number, posify, Pow, Symbol, sympify
from sympy.parsing.sympy_parser import parse_expr

from dimensionful.dimensions import *

# Dictionary holding information of known unit symbols. The key is the symbol,
# the value is a tuple with the conversion factor to cgs, and the
# dimensionality.
unit_symbols_dict = {
    # base
    "g":  (1, mass),
    "cm": (1, length),
    "s":  (1, time),
    "K":  (1, temperature),

    # other cgs
    "dyne": (1, force),
    "erg":  (1, energy),
    "esu":  (1, charge),

    # some SI
    "m": (1e2, length),
    "J": (1e7, energy),
    "Hz": (1, rate),

    # times
    "min": (60, time),
    "hr":  (3600, time),
    "day": (86400, time),     # check cf
    "yr":  (31536000, time),  # check cf

    # Solar units
    "Msun": (1.98892e33, mass),
    "Rsun": (6.96e10, length),
    "Lsun": (3.9e33, power),
    "Tsun": (5870, temperature),

    # astro distances
    "AU": (1.49598e13, length),
    "ly": (9.46053e17, length),
    "pc": (3.08568e18, length),

    # other astro
    "H_0": (2.3e-18, rate),  # check cf

    # other energy units
    "eV": (1.6021766e-12, energy),

    # electric stuff
    "gauss": (1, magnetic_field),
}

# This dictionary formatting from magnitude package, credit to Juan Reyero.
unit_prefixes = {
    'Y': 1e24,   # yotta
    'Z': 1e21,   # zetta
    'E': 1e18,   # exa
    'P': 1e15,   # peta
    'T': 1e12,   # tera
    'G': 1e9,    # giga
    'M': 1e6,    # mega
    'k': 1e3,    # kilo
    'd': 1e1,    # deci
    'c': 1e2,    # centi
    'm': 1e-3,   # mili
    'u': 1e-6,   # micro
    'n': 1e-9,   # nano
    'p': 1e-12,  # pico
    'f': 1e-15,  # femto
    'a': 1e-18,  # atto
    'z': 1e-21,  # zepto
    'y': 1e-24,  # yocto
}

class Unit(Expr):
    """
    Using sympy to represent units as symbols. We just supply extra methods
    so sympy understands how units are related.

    """

    is_positive = True    # make sqrt(m**2) --> m
    is_commutative = True
    is_number = False

    __slots__ = ["expr", "cgs_value", "dimensions", "is_atomic"]

    def __new__(cls, unit_expr=None, cgs_value=None, dimensions=None,
                **assumptions):
        """
        Build a new unit. May be an atomic unit (like a gram) or a combination
        of other units (like g / cm**3). Either way, you can make the unit
        symbol anything.

        Parameters
        ----------
        unit_expr : string or sympy.core.expr.Expr
            The symbolic expression. Symbol("g") for gram.
        cgs_value : float
            This unit's value in cgs. 1.0 for gram.
        dimensions : sympy.core.expr.Expr
            A sympy expression representing the dimensionality of this unit.
            Should just be a sympy.core.mul.Mul object of mass, length, time,
            and temperature objects to various powers. mass for gram.

        """
        # Check for no args
        if not unit_expr:
            unit_expr = sympify(1)

        # if we have a string, parse into an expression
        if isinstance(unit_expr, str):
            unit_expr = parse_expr(unit_expr)

        if not isinstance(unit_expr, Expr):
            raise Exception("Unit representation must be a string or sympy Expr. %s is a %s" % (unit_expr, type(unit_expr)))
        # done with argument checking...

        # sympify, posify, and nsimplify the expr
        unit_expr = sympify(unit_expr)
        p, r = posify(unit_expr)
        unit_expr = p.subs(r)
        unit_expr = nsimplify(unit_expr)

        # see if the unit is atomic.
        is_atomic = False
        if isinstance(unit_expr, Symbol):
            is_atomic = True

        # did they supply cgs_value and dimensions?
        if cgs_value and not dimensions or dimensions and not cgs_value:
            raise Exception("If you provide cgs_vale or dimensions, you must provide both! cgs_value is %s, dimensions is %s." % (cgs_value, dimensions))

        if cgs_value and dimensions:
            # check that cgs_vale is a float or can be converted to one
            try:
                cgs_value = float(cgs_value)
            except ValueError:
                raise ValueError("Please provide a float for the cgs_value kwarg. I got a '%s'." % cgs_value)
            # check that dimensions is valid
            dimensions = verify_dimensions(sympify(dimensions))
            # save the values
            this_cgs_value, this_dimensions = cgs_value, dimensions

        else:  # lookup the unit symbols
            this_cgs_value, this_dimensions = \
                get_unit_data_from_expr(unit_expr)

        # cool trick to get dimensions powers as Rationals
        this_dimensions = nsimplify(this_dimensions)

        # init obj with superclass construct
        obj = Expr.__new__(cls, **assumptions)

        # attach attributes to obj
        obj.expr = unit_expr
        obj.is_atomic = is_atomic
        obj.cgs_value = this_cgs_value
        obj.dimensions = this_dimensions

        # return `obj` so __init__ can handle it.
        return obj

    ### some sympy conventions I guess
    def __getnewargs__(self):
        return (self.expr, self.is_atomic, self.cgs_value, self.dimensions)

    def __hash__(self):
        return super(Unit, self).__hash__()

    def _hashable_content(self):
        return (self.expr, self.is_atomic, self.cgs_value, self.dimensions)
    ### end sympy conventions

    def __repr__(self):
        if self.expr == 1:
            return "(dimensionless)"
        return str(self.expr)

    def __str__(self):
        if self.expr == 1:
            return "(dimensionless)"
        return str(self.expr)

    # for sympy.printing
    def _sympystr(self, *args):
        return str(self.expr)

    ### override sympy operations
    def __mul__(self, right_object):
        """ Multiply Unit with right_object (Unit). """
        return Unit(self.expr * right_object.expr,
                    self.cgs_value * right_object.cgs_value,
                    self.dimensions * right_object.dimensions)

    def __div__(self, right_object):
        """ Divide Unit by right_object (Unit). """
        return Unit(self.expr / right_object.expr,
                    self.cgs_value / right_object.cgs_value,
                    self.dimensions / right_object.dimensions)

    def __pow__(self, power):
        """ Take Unit to a power. """
        return Unit(self.expr**power, self.cgs_value**power,
                    self.dimensions**power)

    ### Comparison operators
    def same_dimensions_as(self, other_unit):
        """ Test if dimensions are the same. """
        return (self.dimensions / other_unit.dimensions) == 1

    def __eq__(self, right_object):
        """
        Test equality. If another Unit, check dimensions and cgs_value. If a
        sympy object, all bets are off.

        """
        if isinstance(right_object, Unit):
            return ( self.cgs_value == right_object.cgs_value
                     and self.dimensions == right_object.dimensions )
        return False

    @property
    def is_dimensionless(self):
        return self.dimensions == 1

    def get_cgs_equivalent(self):
        """ Create and return dimensionally-equivalent cgs units. """
        cgs_units_string = "g**(%s) * cm**(%s) * s**(%s) * K**(%s)" % \
            (self.dimensions.as_coeff_exponent(mass)[1],
             self.dimensions.as_coeff_exponent(length)[1],
             self.dimensions.as_coeff_exponent(time)[1],
             self.dimensions.as_coeff_exponent(temperature)[1])
        return Unit(cgs_units_string, 1, self.dimensions)

# @todo: simpler method that doesn't use recursion would be better...
def verify_dimensions(dimensions):
    """
    Make sure that dimensions is a valid dimension expression. It has to be made
    of only the base dimension symbols, to powers, multiplied together. If
    valid, return the simplified expression. If not, raise an Exception.

    """
    # in the case of a Number of Symbol, we can just return
    if isinstance(dimensions, Number):
        return dimensions
    elif isinstance(dimensions, Symbol):
        if dimensions in base_dimensions:
            return dimensions
        else:
            raise Exception("Dimensions expression contains a non-base dimension symbol '%s'" % str(dimensions))

    # validate args of a Pow or Mul separately
    elif isinstance(dimensions, Pow):
        return verify_dimensions(dimensions.args[0])**verify_dimensions(dimensions.args[1])

    elif isinstance(dimensions, Mul):
        total_mul = 1
        for arg in dimensions.args:
            total_mul *= verify_dimensions(arg)
        return total_mul

    print dimensions

    # should never get here
    raise Exception("Bad dimensions expression.")

def get_unit_data_from_expr(unit_expr):
    """
    Gets total cgs_value and dimensions from a unit expression.

    """
    # sometimes a unit object slips in
    if isinstance(unit_expr, Unit):
        return (unit_expr.cgs_value, unit_expr.dimensions)
    # now for the sympy possibilities
    elif isinstance(unit_expr, Symbol):
        return lookup_unit_symbol(str(unit_expr))

    elif isinstance(unit_expr, Number):
        return (1, 1)

    elif isinstance(unit_expr, Pow):
        unit_data = get_unit_data_from_expr(unit_expr.args[0])
        power = unit_expr.args[1]
        return (unit_data[0]**power, unit_data[1]**power)

    elif isinstance(unit_expr, Mul):
        cgs_value = 1
        dimensions = 1
        for i, expr in enumerate(unit_expr.args):
            unit_data = get_unit_data_from_expr(expr)
            cgs_value *= unit_data[0]
            dimensions *= unit_data[1]

        return (cgs_value, dimensions)

    raise Exception("Cannot parse for unit data from '%s'. Please supply an expression of only Unit/Symbol, Pow, and Mul." % str(unit_expr))

def lookup_unit_symbol(symbol_string):
    """
    Find the unit data of this symbol. Raise an exception if not found.

    """

    if symbol_string in unit_symbols_dict:
        # lookup successful, return the tuple directly
        return unit_symbols_dict[symbol_string]

    # could still be a known symbol with a prefix
    possible_prefix = symbol_string[0]
    if possible_prefix in unit_prefixes:
        # the first character could be a prefix, check the rest of
        # the symbol
        symbol_wo_prefix = symbol_string[1:]

        if symbol_wo_prefix in unit_symbols_dict:
            # lookup successful, it's a symbol with a prefix
            unit_data = unit_symbols_dict[symbol_wo_prefix]
            prefix_value = unit_prefixes[possible_prefix]

            # don't forget to account for the prefix value!
            return (unit_data[0] * prefix_value, unit_data[1])

    # no dice
    raise Exception("Lookup failed. Unknown unit symbol '%s'. Please supply the dimensions and cgs value when creating this object." % symbol_string)

# util function
def get_conversion_factor(old_units, new_units):
    """
    Use the conversion factors table to figure out the factor between these two
    units.

    Parameters
    ----------
    old_units: Unit object
        The current units.
    new_units : Unit object
        The units we want.

    Returns
    -------
    conversion_factor : float
        ``old_units / new_units``
    """
    # @todo: avoid lossy cast
    return float(old_units.cgs_value / new_units.cgs_value)

max_stars_repo_name: caseywstark/dimensionful

max_issues_repo_name: caseywstark/dimensionful

max_forks_repo_name: caseywstark/dimensionful

reason: from sympy


text:
try:
    from . import generic as g
except BaseException:
    import generic as g

TEST_DIM = (10000, 3)

class CacheTest(g.unittest.TestCase):

    def test_hash(self):
        setup = 'import numpy, trimesh;'
        setup += 'd = numpy.random.random((10000,3));'
        setup += 't = trimesh.caching.tracked_array(d)'

        count = 10000

        g.timeit.timeit(setup=setup,
                        stmt='t._modified_m=True;t.md5()',
                        number=count)
        g.timeit.timeit(setup=setup,
                        stmt='t._modified_c=True;t.crc()',
                        number=count)
        g.timeit.timeit(setup=setup,
                        stmt='t._modified_x=True;t.fast_hash()',
                        number=count)

        m = g.get_mesh('featuretype.STL')
        # log result values
        g.log.info('\nResult\nMD5:\n{}\nCRC:\n{}\nXX:\n{}'.format(
            m.vertices.md5(),
            m.vertices.crc(),
            m.vertices.fast_hash()))

    def test_track(self):
        """
        Check to make sure our fancy caching system only changes
        hashes when data actually changes.
        """
        # generate test data and perform numpy operations
        a = g.trimesh.caching.tracked_array(
            g.np.random.random(TEST_DIM))
        modified = []

        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])
        a[0][0] = 10
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        a[0][0] += 0.1
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        a[:10] = g.np.fliplr(a[:10])
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        a[1] = 5
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])
        a[2:] = 2
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        # these operations altered data and
        # the hash SHOULD have changed
        modified = g.np.array(modified)
        assert (g.np.diff(modified, axis=0) != 0).all()

        # now do slice operations which don't alter data
        modified = []
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])
        b = a[[0, 1, 2]]  # NOQA
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])
        c = a[1:]  # NOQA
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])
        # double slice brah
        a = a[::-1][::-1]
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        # these operations should have been cosmetic and
        # the hash should NOT have changed
        modified = g.np.array(modified)
        assert (g.np.diff(modified, axis=0) == 0).all()

        # now change stuff and see if checksums change
        a = g.trimesh.caching.tracked_array([0, 0, 4])
        modified = []

        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])
        a += 10
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        # assign some new data
        a = g.trimesh.caching.tracked_array([.125, 115.32444, 4],
                                            dtype=g.np.float64)

        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        a += [10, 0, 0]
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        a *= 10
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        # itruediv rather than idiv
        a /= 2.0
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        # idiv
        a /= 2.123
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        a -= 1.0
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        # in place floor division :|
        a //= 2
        modified.append([int(a.md5(), 16),
                         a.crc(),
                         a.fast_hash()])

        # these operations altered data and
        # the hash SHOULD have changed
        modified = g.np.array(modified)
        assert (g.np.diff(modified, axis=0) != 0).all()

    def test_contiguous(self):
        a = g.np.random.random((100, 3))
        t = g.trimesh.caching.tracked_array(a)

        # hashing will fail on non- contiguous arrays
        # make sure our utility function has handled this properly
        # for both MD5 and CRC
        assert t.md5() != t[::-1].md5()
        assert t.crc() != t[::-1].crc()
        assert t.fast_hash() != t[::-1].fast_hash()

    def test_mutable(self):
        """
        Run some simple tests on mutable DataStore objects
        """
        d = g.trimesh.caching.DataStore()

        d['hi'] = g.np.random.random(100)
        hash_initial = d.fast_hash()
        # mutate internal data
        d['hi'][0] += 1
        assert d.fast_hash() != hash_initial

        # should be mutable by default
        assert d.mutable
        # set data to immutable
        d.mutable = False

        try:
            d['hi'][1] += 1
        except ValueError:
            # should be raised when array is marked as read only
            return
        # we shouldn't make it past the try-except
        raise ValueError('mutating data worked when it shouldn\'t!')

    def test_transform(self):
        """
        apply_transform tries to not dump the full cache
        """
        m = g.get_mesh('featuretype.STL')
        # should populate edges_face
        e_len = len(m.edges)
        # should maintain required properties
        m.apply_transform(g.transforms[1])
        # should still be in the cache
        assert len(m.edges_face) == e_len

    def test_simple_collision(self):
        faces1 = g.np.array([0, 1, 2, 0, 3, 1, 0,
                             2, 4, 0, 4, 5, 5, 6,
                             3, 5, 3, 0, 7, 1, 3,
                             7, 3, 6, 4, 2, 1, 4,
                             1, 7, 5, 4, 7, 6, 5, 7],
                            dtype=g.np.int64).reshape(-1, 3)
        faces2 = g.np.array([0, 1, 2, 0, 3, 1, 2,
                             4, 0, 5, 4, 2, 6, 3,
                             0, 6, 0, 4, 6, 1, 3,
                             6, 7, 1, 2, 7, 5, 2,
                             1, 7, 4, 5, 7, 6, 4, 7],
                            dtype=g.np.int64).reshape(-1, 3)
        fast_hash = g.trimesh.caching.fast_hash
        assert fast_hash(faces1) != fast_hash(faces2)

if __name__ == '__main__':
    g.trimesh.util.attach_to_log()
    g.unittest.main()

max_stars_repo_name: friggog/trimesh

max_issues_repo_name: friggog/trimesh

max_forks_repo_name: friggog/trimesh


text:
import matplotlib.pyplot as plt
import seaborn as sns
from matplotlib.colors import ListedColormap
from scipy.stats import gaussian_kde
from scipy.ndimage import gaussian_filter
import scipy
import itertools
import numpy as np
import matplotlib as mpl
from scipy.ndimage.filters import gaussian_filter
from .base_plotting import PlotStructure

gray4 = (189/255., 189/255., 189/255.)
gray5 = (150/255., 150/255., 150/255.)
blue2 = (222/255., 235/255., 247/255.)
blue5 = (107/255., 174/255., 214/255.)
orange3 = (253 / 255., 208 / 255., 162 / 255.)
orange5 = (253 / 255., 141 / 255., 60 / 255.)
red5 = (251/255., 106/255., 74/255.)
red6 = (239/255., 59/255., 44/255.)
purple5 = (158/255., 154/255., 200/255.)
purple6 = (128/255., 125/255., 186/255.)
purple9 = (63/255., 0/255., 125/255.)

custom_cmap = ListedColormap([gray4, gray5, blue2, blue5, orange3, orange5,  red5, red6, purple5, purple6, purple9])

class PlotScatter(PlotStructure):
    """
    PlotScatter handles plotting 2D scatter between a set of features.
    It will also optionally overlay KDE contours of the target variable,
    which is a first-order method for evaluating possible feature interactions
    and whether the learned relationships are consistent with the data.
    """
    oranges = ListedColormap(["xkcd:peach", "xkcd:orange", "xkcd:bright orange", "xkcd:rust brown"])
    blues = ListedColormap(["xkcd:periwinkle blue", "xkcd:clear blue", "xkcd:navy blue"])

    def __init__(self, BASE_FONT_SIZE=12):
        super().__init__(BASE_FONT_SIZE=BASE_FONT_SIZE)

    def plot_scatter(self, estimators, X, y, features, display_feature_names={}, display_units={},
                subsample=1.0, peak_val=None, kde=False, **kwargs):
        """
        Plot KDE between two features and color-code by the target variable
        """

        #TODO: Plot relationships for multiple features!!

        estimator_names = list(estimators.keys())
        n_panels = len(estimator_names)
        only_one_model = len(estimator_names)==1

        predictions=np.zeros((n_panels, len(y)), dtype=np.float16)
        j=0
        for estimator_name, estimator in estimators.items():
            if hasattr(estimator, 'predict_proba'):
                predictions[j,:] = estimator.predict_proba(X)[:,1]
            else:
                predictions[j,:] = estimator.predict(X)
            j+=1

        kwargs = self.get_fig_props(n_panels, **kwargs)
        # create subplots, one for each feature
        fig, axes = self.create_subplots(
            n_panels=n_panels,
            sharex=False,
            sharey=False,
            **kwargs,
        )

        ax_iterator = self.axes_to_iterator(n_panels, axes)
        vmax = np.max(predictions)

        for i, ax in enumerate(ax_iterator):
            cf = self.scatter_(ax=ax,
                 X=X,
                 features=features,
                 predictions=predictions[i,:],
                 vmax = vmax,
                 **kwargs,
                         )

            if subsample < 1.0:
                size = min(int(subsample*len(y)), len(y))
                idxs = np.random.choice(len(y), size=size, replace=False)
                var1 = X[features[0]].values[idxs]
                var2 = X[features[1]].values[idxs]
                y1 = y.values[idxs]
            else:
                var1 = X[features[0]].values
                var2 = X[features[1]].values
                y1 = np.copy(y)

            #Shuffle values
            index = np.arange(len(var1))
            np.random.shuffle(index)
            var1 = var1[index]
            var2 = var2[index]
            y1 = y1[index]

            ax.set_xlabel(display_feature_names.get(features[0], features[0]))
            ax.set_ylabel(display_feature_names.get(features[1], features[1]))

            ax.grid(color='#2A3459', alpha=0.6, linestyle='dashed', linewidth=0.5)
            # bluish dark grey, but slightly lighter than background

            if kde:
                cmap_set = [self.oranges, self.blues, 'Reds', 'jet',]
                classes = np.unique(y1)
                idx_sets = [np.where(y1==c) for c in classes]
                for idxs,cmap in zip(idx_sets, cmap_set):
                    # Plot positive cases
                    cs = self.plot_kde_contours(ax,dy=var2[idxs],dx=var1[idxs],
                    target=y1[idxs], cmap=cmap, )
                    handles_set = [cs.legend_elements()[-1]]
                labels = classes

                legend = ax.legend(handles_set, labels, framealpha=0.5)

            # Hide the right and top spines
            ax.spines['right'].set_visible(False)
            ax.spines['top'].set_visible(False)

            if not only_one_model:
                ax.set_title(estimator_names[i])

        if n_panels == 1:
            ax_=axes
        else:
            ax_=axes.ravel().tolist()

        if hasattr(estimator,'predict_proba'):
            cbar_label = 'Probability'
        else:
            cbar_label = 'Response'

        cbar_label = kwargs.get('cbar_label', cbar_label)

        fig.colorbar(cf,ax=ax_,label=cbar_label, orientation='horizontal')

        return fig, axes

    def scatter_(self, ax, features, X, predictions, vmax, **kwargs):
        """
        Plot 2D scatter of ML predictions;
        Only plots a random 20000 points
        """
        size = min(20000, len(X))
        idxs = np.random.choice(len(X), size=size, replace=False)

        x_val = X[features[0]].values[idxs]
        y_val = X[features[1]].values[idxs]
        z_val = predictions[idxs]

        # Show highest predictions on top!
        index = np.argsort(z_val)
        #index = np.random.choice(len(z_val), size=len(z_val), replace=False)

        x_val = x_val[index]
        y_val = y_val[index]
        z_val = z_val[index]

        cmap = kwargs.get('cmap' , custom_cmap)
        zmax = vmax+0.05 if vmax < 1.0 else 1.1
        delta = 0.05 if vmax < 0.5 else 0.1
        levels= [0, 0.05] + list(np.arange(0.1, zmax, delta))
        norm = mpl.colors.BoundaryNorm(levels, cmap.N)

        sca = ax.scatter(x_val,
                     y_val,
                     c=z_val,
                     cmap=cmap,
                     alpha=0.6,
                     s=3,
                     norm=norm,
                    )
        return sca

    def kernal_density_estimate(self, dy, dx):
        dy_min = np.amin(dy)
        dx_min = np.amin(dx)
        dy_max = np.amax(dy)
        dx_max = np.amax(dx)

        x, y = np.mgrid[dx_min:dx_max:100j, dy_min:dy_max:100j]
        positions = np.vstack([x.ravel(), y.ravel()])
        values = np.vstack([dx, dy])
        kernel = gaussian_kde(values)
        f = np.reshape(kernel(positions).T, x.shape)

        return x, y, f

    def plot_kde_contours(self, ax,dy,dx,target,cmap,):
        x, y, f = self.kernal_density_estimate(dy,dx)
        temp_linewidths = [ 0.85, 1.0, 1.25, 1.75]
        temp_thresh = [75.0, 90., 95., 97.5]
        temp_levels = [0.0, 0.0, 0.0, 0.0]
        for i in range(0, len(temp_thresh)):
            temp_levels[i] = np.percentile(f.ravel(), temp_thresh[i])

        #masked_f = np.ma.masked_where(f < 1.6e-5, f)
        cs = ax.contour(
            x,
            y,
            f,
            levels=temp_levels,
            cmap = cmap,
            linewidths=temp_linewidths,
            alpha=1.0,
        )
        fmt = {}
        for l, s in zip(cs.levels, temp_thresh[::-1]):
            fmt[l] = f'{int(s)}%'
        ax.clabel(cs, cs.levels, inline=True, fontsize=self.FONT_SIZES['teensie'], fmt=fmt)

        return cs

max_stars_repo_name: monte-flora/py-mint

max_issues_repo_name: monte-flora/py-mint

max_forks_repo_name: monte-flora/mintpy

reason: import numpy,import scipy,from scipy