text:
#! /usr/bin/python
# -*- coding: utf-8 -*-

# import funkcí z jiného adresáře
import sys
import os.path
import io3d
path_to_script = os.path.dirname(os.path.abspath(__file__))
#import featurevector
import unittest

from loguru import logger
# logger = logging.getLogger()

# import apdb
# apdb.set_trace();
# import scipy.io
import numpy as np

# ----------------- my scripts --------
import argparse
from . import segmentation_general

class VesselSegmentation:
    def __init__(self):
        pass

    def set_params(self):
        pass

    def run(self):
        pass

    def get_output(self):
        pass

if __name__ == "__main__":

    ## logger = logging.getLogger()
    # logger = logging.getLogger()

    logger.setLevel(logging.WARNING)
    ch = logging.StreamHandler()
    logger.addHandler(ch)

    #logger.debug('input params')

    # input parser
    parser = argparse.ArgumentParser(description='Segment vessels from liver')
    parser.add_argument('-dd','--dcmdir',
            default=None,
            help='path to data dir')
    parser.add_argument('-d', '--debug', action='store_true',
            help='run in debug mode')
    parser.add_argument('-i', '--inputfile',  default=None,
            help='input file or directory with data')
    args = parser.parse_args()

    if args.debug:
        logger.setLevel(logging.DEBUG)

    defaultoutputfile =  "vessels.pkl"
    if args.defaultoutputfile:
        args.outputfile = defaultoutputfile

    #else:
    #dcm_read_from_dir('/home/mjirik/data/medical/data_orig/46328096/')
        #data3d, metadata = dcmreaddata.dcm_read_from_dir()
    datap = io3d.read(args.inputfile)

    import sed3
    # pyed = sed3.sed3(oseg.orig_scale_segmentation)
    #pyed.show()
# information about crop
    #cri = oseg.crinfo
    #oseg.data3d = oseg.data3d[cri[0][0]:cri[0][1],cri[1][0]:cri[1][1],cri[2][0]:cri[2][1]]
    #pyed = sed3.sed3(oseg.data3d, contour = oseg.orig_scale_segmentation)

    print('slab', datap['slab'])
    #import ipdb; ipdb.set_trace()  # BREAKPOINT
    #pyed = sed3.sed3(data['data3d'], contour = data['segmentation'])
    #pyed.show()
    #import pdb; pdb.set_trace()

    outputTmp = segmentation_general.vesselSegmentation(
        datap['data3d'],
        segmentation = datap['segmentation'],
        #segmentation = oseg.orig_scale_segmentation,
        threshold = -1,
        inputSigma = 0.15,
        dilationIterations = 2,
        nObj = 1,
        biggestObjects = args.biggest,
#        dataFiltering = True,
        interactivity = True,
        binaryClosingIterations = 2,
        binaryOpeningIterations = 0)

    datap['slab']['none'] = 0
    datap['slab']['liver'] = 1
    datap['slab']['porta'] = 2

    #print np.max(output)
    #import pdb; pdb.set_trace()
    #data = {}
    #data['data3d'] = oseg.data3d
    #data['crinfo'] = oseg.crinfo
    #data['segmentation'] = oseg.segmentation
    datap['segmentation'][output] = datap['slab']['porta']
    #data['slab'] = slab

    pyed = sed3.sed3(datap['data3d'], contour=datap['segmentation'] == datap['slab']['porta'])
    pyed.show()

    #pyed = sed3.sed3(data['segmentation'])
    #pyed.show()
    # Uvolneni pameti

    if args.outputfile == None:
        io3d.write(datap, args.outpufile)

max_stars_repo_name: mjirik/lisa

max_issues_repo_name: mjirik/lisa

max_forks_repo_name: mjirik/lisa


text:
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.collections import PolyCollection
from bokeh import mpl

# Generate data. In this case, we'll make a bunch of center-points and generate
# verticies by subtracting random offsets from those center-points
numpoly, numverts = 100, 4
centers = 100 * (np.random.random((numpoly, 2)) - 0.5)
offsets = 10 * (np.random.random((numverts, numpoly, 2)) - 0.5)
verts = centers + offsets
verts = np.swapaxes(verts, 0, 1)

# In your case, "verts" might be something like:
# verts = zip(zip(lon1, lat1), zip(lon2, lat2), ...)
# If "data" in your case is a numpy array, there are cleaner ways to reorder
# things to suit.

facecolors = ['red', 'green', 'blue', 'cyan', 'yellow', 'magenta', 'black']

edgecolors = ['cyan', 'yellow', 'magenta', 'black', 'red', 'green', 'blue']

widths = [5, 10, 20, 10, 5]

ax = plt.axes()

# Make the collection and add it to the plot.
col = PolyCollection(verts, facecolor=facecolors, edgecolor=edgecolors,
                     linewidth=widths, linestyle='--', alpha=0.5)
ax.add_collection(col)

plt.xlim([-60, 60])
plt.ylim([-60, 60])

plt.title("MPL-PolyCollection support in Bokeh")

mpl.to_bokeh(name="polycollection")

ext: py

max_stars_repo_name: tswicegood/bokeh

reason: import numpy


text:
from numpy import arange, zeros
from codetools.contexts.with_mask import Mask

dep_len = depth.shape
vp = zeros(dep_len)
vs = zeros(dep_len)
with Mask((depth < 4000.0) & (depth > 1000.0)):
    vp = arange(0., 10., 1.)
    vs = arange(0., 100., 10.)
## mgr = Mask((depth < 4000.0) & (depth > 1000.0))
## exit_code = mgr.__exit__
## mgr.__enter__()
## exc = True
## try:
##     try:
##         vp = arange(0., 10., 1.)
##         vs = arange(0., 100., 10.)
##     except:
##         exc = False
## finally:
##     if exc:
##         exit_code(None, None, None)

reason: from numpy


text:
import itertools
from typing import Dict, Tuple, Iterable

import numpy as np
import xarray as xr

from xcube.util.dsio import FORMAT_NAME_ZARR, FORMAT_NAME_NETCDF4

def chunk_dataset(dataset: xr.Dataset,
                  chunk_sizes: Dict[str, int] = None,
                  format_name: str = None) -> xr.Dataset:
    """
    Chunk dataset and update encodings for given format.

    :param dataset: input dataset
    :param chunk_sizes: mapping from dimension name to new chunk size
    :param format_name: format, e.g. "zarr" or "netcdf4"
    :return: the re-chunked dataset
    """
    chunked_ds = dataset.chunk(chunks=chunk_sizes)

    # Update encoding so writing of chunked_ds recognizes new chunks
    chunk_sizes_attr_name = None
    if format_name == FORMAT_NAME_ZARR:
        chunk_sizes_attr_name = "chunks"
    if format_name == FORMAT_NAME_NETCDF4:
        chunk_sizes_attr_name = "chunksizes"
    if chunk_sizes_attr_name:
        for var_name in chunked_ds.variables:
            var = chunked_ds[var_name]
            if chunk_sizes:
                sizes = tuple(chunk_sizes[dim_name] if dim_name in chunk_sizes
                              else var.shape[var.dims.index(dim_name)]
                              for dim_name in var.dims)
                var.encoding.update({chunk_sizes_attr_name: sizes})
            elif chunk_sizes_attr_name in var.encoding:
                # Remove any explicit and wrong specification so writing will use Dask chunks (TBC!)
                del var.encoding[chunk_sizes_attr_name]

    return chunked_ds

def get_empty_dataset_chunks(dataset: xr.Dataset) -> Dict[str, Tuple[Tuple[int, ...]]]:
    """
    Identify empty dataset chunks and return their indices.

    :param dataset: The dataset.
    :return: A mapping from variable name to a list of block indices.
    """
    return {var_name: get_empty_var_chunks(dataset[var_name]) for var_name in dataset.data_vars}

def get_empty_var_chunks(var: xr.DataArray) -> Tuple[Tuple[int, ...]]:
    """
    Identify empty variable chunks and return their indices.

    :param var: The variable.
    :return: A list of block indices.
    """
    chunks = var.chunks
    if chunks is None:
        raise ValueError('data array not chunked')

    chunk_slices = compute_chunk_slices(chunks)

    empty_chunk_indexes = []
    for chunk_index, chunk_slice in chunk_slices:
        data_index = tuple(slice(start, end) for start, end in chunk_slice)
        data = var[data_index]
        if np.all(np.isnan(data)):
            empty_chunk_indexes.append(chunk_index)
            # print(f'empty: {var.name}/{".".join(map(str, chunk_index))}')

    # noinspection PyTypeChecker
    return tuple(empty_chunk_indexes)

def compute_chunk_slices(chunks: Tuple[Tuple[int, ...], ...]) -> Iterable:

    chunk_indices = []
    for c in chunks:
        chunk_indices.append(tuple(i for i in range(len(c))))

    chunk_slices = []
    for c in chunks:
        x = []
        o = 0
        for s in c:
            x.append((o, o + s))
            o += s
        chunk_slices.append(tuple(x))
    return zip(itertools.product(*chunk_indices), itertools.product(*chunk_slices))

max_stars_repo_name: dzelge/xcube

max_issues_repo_name: dzelge/xcube

max_forks_repo_name: dzelge/xcube


text:
#!/usr/bin/env python3

from __future__ import division, print_function, absolute_import

import roslib
import rospy
import numpy as np
import cv2
from qcar.q_essential import LIDAR

from sensor_msgs.msg import LaserScan
from std_msgs.msg import String
import time

class LIDARNode(object):
	def __init__(self):
		super().__init__()
		self.lidar_pub = rospy.Publisher('/scan', LaserScan, queue_size=1000)
		self.num_measurements = 720
		self.lidar = LIDAR(num_measurements=self.num_measurements)
		self.distances = np.zeros((self.num_measurements,1))
		self.angles = np.zeros((self.num_measurements,1))
		self.sampleTime = 1 / 30
		#Get readings
		while not rospy.is_shutdown():
			# sampleTime = 1 / 30
			starTime = time.time()
			self.lidar.read()
			mid_time = time.time()
			scan_time = mid_time - starTime
			# print(scan_time)
			# Calculate the computation time, and the time that the thread should pause/sleep for
			computationTime = scan_time
			sleepTime = self.sampleTime - ( computationTime % self.sampleTime )

			# Pause/sleep and print out the current timestamp
			time.sleep(sleepTime)
			# counter += 1
			self.process_lidar_data(self.lidar_pub, self.lidar.distances.astype(np.float32), self.num_measurements, scan_time)
			process_time = time.time() - mid_time
			# print(process_time, scan_time)
		self.lidar.terminate()

#--------------------------------------------------------------------------------------------------------------
	def process_lidar_data(self, pub_lidar, distances, num_measurement, scan_times):

		scan = LaserScan()
		scan.header.stamp = rospy.Time.now()
		scan.header.frame_id = 'lidar'
		scan.angle_min = 0.0
		scan.angle_max = 6.2744
		# scan.angle_increment = 6.2744 / num_measurement
		# scan.time_increment = scan_times / num_measurement
		scan.angle_increment = 0.0174532923847
		scan.time_increment = 0.000132342218421
		scan.scan_time = scan_times
		scan.range_min = 0.15
		scan.range_max = 12.0
		# print('before')
		scan.ranges = distances.tolist()
		# scan.ranges = []
		# print('after')
		# for i in range(num_measurement):
		# 	scan.ranges.append(distances[i])
		pub_lidar.publish(scan)

if __name__ == '__main__':
	rospy.init_node('lidar_node')
	r = LIDARNode()

	rospy.spin()

max_stars_repo_name: bchampp/scylla

max_issues_repo_name: bchampp/scylla

max_forks_repo_name: bchampp/scylla


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

import numpy

def cg_spin(l, jdouble, mjdouble, spin):
    '''Clebsch Gordon coefficient of <l,m,1/2,spin|j,mj>'''
    ll1 = 2 * l + 1
    if jdouble == 2*l+1:
        if spin > 0:
            c = numpy.sqrt(.5*(ll1+mjdouble)/ll1)
        else:
            c = numpy.sqrt(.5*(ll1-mjdouble)/ll1)
    elif jdouble == 2*l-1:
        if spin > 0:
            c =-numpy.sqrt(.5*(ll1-mjdouble)/ll1)
        else:
            c = numpy.sqrt(.5*(ll1+mjdouble)/ll1)
    else:
        c = 0
    return c

if __name__ == '__main__':
    for kappa in list(range(-4,0)) + list(range(1,4)):
        if kappa < 0:
            l = -kappa - 1
            j = l * 2 + 1
        else:
            l = kappa
            j = l * 2 - 1
        print(kappa,l,j)
        for mj in range(-j, j+1, 2):
            print(cg_spin(l, j, mj, 1), cg_spin(l, j, mj, -1))

max_stars_repo_name: robert-anderson/pyscf

max_issues_repo_name: robert-anderson/pyscf

max_forks_repo_name: robert-anderson/pyscf


text:
import argparse
import os

from osgeo import gdal
import h5py
import numpy as np

import isce3
from nisar.workflows import geocode_insar, h5_prep
from nisar.workflows.geocode_insar_runconfig import GeocodeInsarRunConfig

import iscetest

def test_geocode_run():
    '''
    Test if phase geocodeping runs smoothly
    '''
    test_yaml_path = os.path.join(iscetest.data, 'insar_test.yaml')
    for pu in ['cpu', 'gpu']:
        # Skip GPU geocode insar if cuda not included
        if pu == 'gpu' and not hasattr(isce3, 'cuda'):
            continue

        with open(test_yaml_path) as fh_test_yaml:
            test_yaml = fh_test_yaml.read().replace('@ISCETEST@', iscetest.data). \
                replace('@TEST_OUTPUT@', f'{pu}_gunw.h5'). \
                replace('@TEST_PRODUCT_TYPES@', 'GUNW'). \
                replace('@TEST_RDR2GEO_FLAGS@', 'True')
            if pu == 'gpu':
                test_yaml = test_yaml.replace('gpu_enabled: False', 'gpu_enabled: True')

        # Create CLI input name space with yaml instead of filepath
        args = argparse.Namespace(run_config_path=test_yaml, log_file=False)

        # Initialize runconfig object
        runconfig = GeocodeInsarRunConfig(args)
        runconfig.geocode_common_arg_load()

        # save geogrid for later validation
        geogrid = runconfig.cfg['processing']['geocode']['geogrids']['A']
        geo_tx = np.array([geogrid.start_x, geogrid.spacing_x, 0, geogrid.start_y, 0, geogrid.spacing_y])
        geo_tx.tofile('gunw_geogrid.txt', sep=',')

        # prepare HDF5 outputs
        out_paths = h5_prep.run(runconfig.cfg)

        # insert rdr2geo outputs into RUNW HDF5
        rdr2geo_dict = {'x': 'unwrappedPhase', 'y': 'coherenceMagnitude'}
        runw_product_path = 'science/LSAR/RUNW/swaths/frequencyA/interferogram/HH'
        az_looks = runconfig.cfg['processing']['crossmul']['azimuth_looks']
        rg_looks = runconfig.cfg['processing']['crossmul']['range_looks']
        with h5py.File(out_paths['RUNW'], 'a', libver='latest', swmr=True) as h_runw:
            for axis, ds_name in rdr2geo_dict.items():
                # extract rdr2geo raster to array
                f_rdr2geo = f'rdr2geo/freqA/{axis}.rdr'
                ds = gdal.Open(f_rdr2geo, gdal.GA_ReadOnly)
                arr = ds.GetRasterBand(1).ReadAsArray()
                ds = None

                # decimate array (take center pixel of mlook grid) and save to HDF5
                sz_x, sz_y = arr.shape
                arr_mlook = arr[az_looks//2:sz_x-az_looks//2:az_looks, az_looks//2:sz_y-az_looks//2:az_looks]
                h_runw[f'{runw_product_path}/{ds_name}'][:, :] = arr_mlook

        # disable unused runw dataset
        runconfig.cfg['processing']['geocode']['datasets']['connected_components'] = False

        # run geocodeing of runw
        geocode_insar.run(runconfig.cfg, out_paths['RUNW'], out_paths['GUNW'])

def test_geocode_validate():
    '''
    validate generated geocode data
    '''
    # load geotransform
    geo_trans = np.fromfile('gunw_geogrid.txt', sep=',')

    # prepare the check gunw hdf5
    scratch_path = '.'
    for pu in ['cpu', 'gpu']:
        # Skip GPU geocode insar if cuda not included
        if pu == 'gpu' and not hasattr(isce3, 'cuda'):
            continue

        path_gunw = os.path.join(scratch_path, f'{pu}_gunw.h5')
        product_path = 'science/LSAR/GUNW/grids/frequencyA/interferogram/HH'
        with h5py.File(path_gunw, 'r', libver='latest', swmr=True) as h:
            # iterate over axis
            rdr2geo_dict = {'x': 'unwrappedPhase', 'y': 'coherenceMagnitude'}
            for axis, ds_name in rdr2geo_dict.items():
                # get dataset as masked array
                ds_path = f'{product_path}/{ds_name}'
                geo_arr = h[ds_path][()]
                geo_arr = np.ma.masked_array(geo_arr, mask=np.isnan(geo_arr))

                # get transform and meshgrids once for common geogrid
                if ds_name == 'unwrappedPhase':
                    x0 = geo_trans[0] + geo_trans[1] / 2.0
                    dx = geo_trans[1]
                    y0 = geo_trans[3] + geo_trans[5] / 2.0
                    dy = geo_trans[5]

                    lines, pixels = geo_arr.shape
                    meshx, meshy = np.meshgrid(np.arange(pixels), np.arange(lines))
                    grid_lon = x0 + meshx * dx
                    grid_lat = y0 + meshy * dy

                # set calculate error and RMSE error thresholds per axis
                if ds_name == 'unwrappedPhase':
                    err = geo_arr - grid_lon
                    rmse_err_threshold = 0.5 * dx
                else:
                    err = geo_arr - grid_lat
                    rmse_err_threshold = 0.5 * abs(dy)

                # check max error
                max_err = np.nanmax(err)
                # 1e-4 is ~11m
                assert(max_err < 1e-4), f'{axis}-axis max error fail on {pu}'

                # check RMSE
                rmse = np.sqrt(np.sum(err**2) / np.count_nonzero(~geo_arr.mask))
                assert(rmse < rmse_err_threshold), f'{axis}-axis RMSE fail on {pu}'

if __name__ == "__main__":
    test_geocode_run()


text:
"""Operations for dannce."""
import numpy as np
import cv2
import time
from typing import Text
import torch
import torch.nn.functional as F

class Camera:
    def __init__(self, R, t, K, tdist, rdist, name=""):
        self.R = np.array(R).copy()
        assert self.R.shape == (3, 3)

        self.t = np.array(t).copy()
        assert self.t.shape == (1, 3)

        self.K = np.array(K).copy()
        assert self.K.shape == (3, 3)

        self.M = np.concatenate((R, t), axis=0) @ self.K

        self.tdist = tdist
        self.rdist = rdist

        self.name = name

    def update_after_crop(self, bbox):
        left, upper, right, lower = bbox

        cx, cy = self.K[2, 0], self.K[2, 1]

        new_cx = cx - left
        new_cy = cy - upper

        self.K[2, 0], self.K[2, 1] = new_cx, new_cy

    def update_after_resize(self, image_shape, new_image_shape):
        height, width = image_shape
        new_height, new_width = new_image_shape

        fx, fy, cx, cy = self.K[0, 0], self.K[1, 1], self.K[2, 0], self.K[2, 1]

        new_fx = fx * (new_width / width)
        new_fy = fy * (new_height / height)
        new_cx = cx * (new_width / width)
        new_cy = cy * (new_height / height)

        self.K[0, 0], self.K[1, 1], self.K[2, 0], self.K[2, 1] = new_fx, new_fy, new_cx, new_cy

    @property
    def camera_matrix(self):
        return self.extrinsics.dot(self.K)

    @property
    def extrinsics(self):
        return np.concatenate((self.R, self.t), axis=0)

def camera_matrix(K: np.ndarray, R: np.ndarray, t: np.ndarray) -> np.ndarray:
    """Derive the camera matrix.

    Derive the camera matrix from the camera intrinsic matrix (K),
    and the extrinsic rotation matric (R), and extrinsic
    translation vector (t).

    Note that this uses the matlab convention, such that
    M = [R;t] * K
    """
    return np.concatenate((R, t), axis=0) @ K

def world_to_cam(pts, M, device):
    M = M.to(device=device)
    pts1 = torch.ones(pts.shape[0], 1, dtype=torch.float32, device=device)

    projPts = torch.matmul(torch.cat((pts, pts1), 1), M)
    return projPts

def project_to2d(pts, M: np.ndarray, device: Text) -> torch.Tensor:
    """Project 3d points to 2d.

    Projects a set of 3-D points, pts, into 2-D using the camera intrinsic
    matrix (K), and the extrinsic rotation matric (R), and extrinsic
    translation vector (t). Note that this uses the matlab
    convention, such that
    M = [R;t] * K, and pts2d = pts3d * M
    """

    # pts = torch.Tensor(pts.copy()).to(device)
    M = M.to(device=device)
    pts1 = torch.ones(pts.shape[0], 1, dtype=torch.float32, device=device)

    projPts = torch.matmul(torch.cat((pts, pts1), 1), M)
    projPts[:, :2] = projPts[:, :2] / projPts[:, 2:]

    return projPts

def sample_grid_nearest(
    im: np.ndarray, projPts: np.ndarray, device: Text
) -> torch.Tensor:
    """Unproject features."""
    # im_x, im_y are the x and y coordinates of each projected 3D position.
    # These are concatenated here for every image in each batch,

    feats = torch.as_tensor(im.copy(), device=device)
    grid = projPts
    c = int(round(projPts.shape[0] ** (1 / 3.0)))

    fh, fw, fdim = list(feats.shape)

    # # make sure all projected indices fit onto the feature map
    im_x = torch.clamp(grid[:, 0], 0, fw - 1)
    im_y = torch.clamp(grid[:, 1], 0, fh - 1)

    im_xr = im_x.round().type(torch.long)
    im_yr = im_y.round().type(torch.long)
    im_xr[im_xr < 0] = 0
    im_yr[im_yr < 0] = 0
    Ir = feats[im_yr, im_xr]

    return Ir.reshape((c, c, c, -1)).permute(3, 0, 1, 2).unsqueeze(0)

def sample_grid_linear(
    im: np.ndarray, projPts: np.ndarray, device: Text
) -> torch.Tensor:
    """Unproject features."""
    # im_x, im_y are the x and y coordinates of each projected 3D position.
    # These are concatenated here for every image in each batch,

    feats = torch.as_tensor(im.copy(), device=device)
    grid = projPts
    c = int(round(projPts.shape[0] ** (1 / 3.0)))

    fh, fw, fdim = list(feats.shape)

    # # make sure all projected indices fit onto the feature map
    im_x = torch.clamp(grid[:, 0], 0, fw - 1)
    im_y = torch.clamp(grid[:, 1], 0, fh - 1)

    # round all indices
    im_x0 = torch.floor(im_x).type(torch.long)
    # new array with rounded projected indices + 1
    im_x1 = im_x0 + 1
    im_y0 = torch.floor(im_y).type(torch.long)
    im_y1 = im_y0 + 1

    # Convert from int to float -- but these are still round
    # numbers because of rounding step above
    im_x0_f, im_x1_f = im_x0.type(torch.float), im_x1.type(torch.float)
    im_y0_f, im_y1_f = im_y0.type(torch.float), im_y1.type(torch.float)

    # Gather  values
    # Samples all featuremaps at the projected indices,
    # and their +1 counterparts. Stop at Ia for nearest neighbor interpolation.

    # need to clip the corner indices because they might be out of bounds...
    # This could lead to different behavior compared to TF/numpy, which return 0
    # when an index is out of bounds
    im_x1_safe = torch.clamp(im_x1, 0, fw - 1)
    im_y1_safe = torch.clamp(im_y1, 0, fh - 1)

    im_x1[im_x1 < 0] = 0
    im_y1[im_y1 < 0] = 0
    im_x0[im_x0 < 0] = 0
    im_y0[im_y0 < 0] = 0
    im_x1_safe[im_x1_safe < 0] = 0
    im_y1_safe[im_y1_safe < 0] = 0

    Ia = feats[im_y0, im_x0]
    Ib = feats[im_y0, im_x1_safe]
    Ic = feats[im_y1_safe, im_x0]
    Id = feats[im_y1_safe, im_x1_safe]

    # To recaptiulate behavior  in numpy/TF, zero out values that fall outside bounds
    Ib[im_x1 > fw - 1] = 0
    Ic[im_y1 > fh - 1] = 0
    Id[(im_x1 > fw - 1) | (im_y1 > fh - 1)] = 0
    # Calculate bilinear weights
    # We've now sampled the feature maps at corners around the projected values
    # Here, the corners are weighted by distance from the projected value
    wa = (im_x1_f - im_x) * (im_y1_f - im_y)
    wb = (im_x1_f - im_x) * (im_y - im_y0_f)
    wc = (im_x - im_x0_f) * (im_y1_f - im_y)
    wd = (im_x - im_x0_f) * (im_y - im_y0_f)

    Ibilin = (
        wa.unsqueeze(1) * Ia
        + wb.unsqueeze(1) * Ib
        + wc.unsqueeze(1) * Ic
        + wd.unsqueeze(1) * Id
    )

    return Ibilin.reshape((c, c, c, -1)).permute(3, 0, 1, 2).unsqueeze(0)

def sample_grid(im: np.ndarray, projPts: np.ndarray, device: Text, method: Text = "linear"):
    """Transfer 3d features to 2d by projecting down to 2d grid, using torch.

    Use 2d interpolation to transfer features to 3d points that have
    projected down onto a 2d grid
    Note that function expects proj_grid to be flattened, so results should be
    reshaped after being returned
    """
    if method == "nearest" or method == "out2d":
        proj_rgb = sample_grid_nearest(im, projPts, device)
    elif method == "linear" or method == "bilinear":
        proj_rgb = sample_grid_linear(im, projPts, device)
    else:
        raise Exception("{} not a valid interpolation method".format(method))

    return proj_rgb

def unDistortPoints(
    pts,
    intrinsicMatrix,
    radialDistortion,
    tangentDistortion,
    rotationMatrix,
    translationVector,
):
    """Remove lens distortion from the input points.

    Input is size (M,2), where M is the number of points
    """
    dcoef = radialDistortion.ravel()[:2].tolist() + tangentDistortion.ravel().tolist()

    if len(radialDistortion.ravel()) == 3:
        dcoef = dcoef + [radialDistortion.ravel()[-1]]
    else:
        dcoef = dcoef + [0]

    ts = time.time()
    pts_u = cv2.undistortPoints(
        np.reshape(pts, (-1, 1, 2)).astype("float32"),
        intrinsicMatrix.T,
        np.array(dcoef),
        P=intrinsicMatrix.T,
    )

    pts_u = np.reshape(pts_u, (-1, 2))

    return pts_u

def triangulate(pts1, pts2, cam1, cam2):
    """Return triangulated 3- coordinates.

    Following Matlab convetion, given lists of matching points, and their
    respective camera matrices, returns the triangulated 3- coordinates.
    pts1 and pts2 must be Mx2, where M is the number of points with
    (x,y) positions. M 3-D points will be returned after triangulation
    """
    pts1 = pts1.T
    pts2 = pts2.T

    cam1 = cam1.T
    cam2 = cam2.T

    out_3d = np.zeros((3, pts1.shape[1]))

    for i in range(out_3d.shape[1]):
        if ~np.isnan(pts1[0, i]):
            pt1 = pts1[:, i : i + 1]
            pt2 = pts2[:, i : i + 1]

            A = np.zeros((4, 4))
            A[0:2, :] = pt1 @ cam1[2:3, :] - cam1[0:2, :]
            A[2:, :] = pt2 @ cam2[2:3, :] - cam2[0:2, :]

            u, s, vh = np.linalg.svd(A)
            v = vh.T

            X = v[:, -1]
            X = X / X[-1]

            out_3d[:, i] = X[0:3].T
        else:
            out_3d[:, i] = np.nan

    return out_3d

def triangulate_multi_instance(pts, cams):
    """Return triangulated 3- coordinates.

    Following Matlab convetion, given lists of matching points, and their
    respective camera matrices, returns the triangulated 3- coordinates.
    pts1 and pts2 must be Mx2, where M is the number of points with
    (x,y) positions. M 3-D points will be returned after triangulation
    """
    pts = [pt.T for pt in pts]
    cams = [c.T for c in cams]
    out_3d = np.zeros((3, pts[0].shape[1]))
    # traces = np.zeros((out_3d.shape[1],))

    for i in range(out_3d.shape[1]):
        if ~np.isnan(pts[0][0, i]):
            p = [p[:, i : i + 1] for p in pts]

            A = np.zeros((2 * len(cams), 4))
            for j in range(len(cams)):
                A[(j) * 2 : (j + 1) * 2] = p[j] @ cams[j][2:3, :] - cams[j][0:2, :]

            u, s, vh = np.linalg.svd(A)
            v = vh.T

            X = v[:, -1]
            X = X / X[-1]

            out_3d[:, i] = X[0:3].T
            # traces[i] = np.sum(s[0:3])

        else:
            out_3d[:, i] = np.nan

    return out_3d

def ravel_multi_index(I, J, shape):
    """Create an array of flat indices from coordinate arrays.

    shape is (rows, cols)
    """
    r, c = shape
    return I * c + J

def distortPoints(
    points, intrinsicMatrix, radialDistortion, tangentialDistortion, device
):
    """Distort points according to camera parameters.
    Ported from Matlab 2018a
    """

    # unpack the intrinsic matrix
    cx = intrinsicMatrix[2, 0]
    cy = intrinsicMatrix[2, 1]
    fx = intrinsicMatrix[0, 0]
    fy = intrinsicMatrix[1, 1]
    skew = intrinsicMatrix[1, 0]

    # center the points
    center = torch.as_tensor((cx, cy), dtype=torch.float32, device=device)
    centeredPoints = points - center

    # normalize the pcenteredPoints[:, 1] / fyoints
    yNorm = centeredPoints[:, 1] / fy
    xNorm = (centeredPoints[:, 0] - skew * yNorm) / fx

    # compute radial distortion
    r2 = xNorm ** 2 + yNorm ** 2
    r4 = r2 * r2
    r6 = r2 * r4

    k = np.zeros((3,))
    k[:2] = radialDistortion[:2]
    if list(radialDistortion.shape)[0] < 3:
        k[2] = 0
    else:
        k[2] = radialDistortion[2]
    alpha = k[0] * r2 + k[1] * r4 + k[2] * r6

    # compute tangential distortion
    p = tangentialDistortion
    xyProduct = xNorm * yNorm
    dxTangential = 2 * p[0] * xyProduct + p[1] * (r2 + 2 * xNorm ** 2)
    dyTangential = p[0] * (r2 + 2 * yNorm ** 2) + 2 * p[1] * xyProduct

    # apply the distortion to the points
    normalizedPoints = torch.transpose(torch.stack((xNorm, yNorm)), 0, 1)

    distortedNormalizedPoints = (
        normalizedPoints
        + normalizedPoints * torch.transpose(torch.stack((alpha, alpha)), 0, 1)
        + torch.transpose(torch.stack((dxTangential, dyTangential)), 0, 1)
    )

    distortedPointsX = (
        distortedNormalizedPoints[:, 0] * fx
        + cx
        + skew * distortedNormalizedPoints[:, 1]
    )

    distortedPointsY = distortedNormalizedPoints[:, 1] * fy + cy

    distortedPoints = torch.stack((distortedPointsX, distortedPointsY))

    return distortedPoints

def expected_value_3d(prob_map, grid_centers):
    bs, channels, h, w, d = prob_map.shape

    prob_map = prob_map.permute(0, 2, 3, 4, 1).reshape(-1, channels)
    grid_centers = grid_centers.reshape(-1, 3)
    weighted_centers = prob_map.unsqueeze(1) * grid_centers.unsqueeze(-1)
    weighted_centers = weighted_centers.reshape(-1, h*w*d, 3, channels).sum(1)

    return weighted_centers # [bs, 3, channels]

def spatial_softmax(feats):
    bs, channels, h, w, d = feats.shape
    feats = feats.reshape(bs, channels, -1)
    feats = F.softmax(feats, dim=-1)
    return feats.reshape(bs, channels, h, w, d)

def var_3d(prob_map, grid_centers, markerlocs):
    """Return the average variance across all marker probability maps.

    Used a loss to promote "peakiness" in the probability map output
    prob_map should be (batch_size,h,w,d,channels)
    grid_centers should be (batch_size,h*w*d,3)
    markerlocs is (batch_size,3,channels)
    """
    channels, h, w, d = prob_map.shape[1:]
    prob_map = prob_map.permute(0, 2, 3, 3, 1).reshape(-1, channels)
    grid_dist = (grid_centers.unsqueeze(-1) - markerlocs.unsqueeze(1)) ** 2
    grid_dist = grid_dist.sum(2)
    grid_dist = grid_dist.reshape(-1, channels)

    weighted_var = prob_map * grid_dist
    weighted_var = weighted_var.reshape(-1, h*w*d, channels)
    weighted_var = weighted_var.sum(1)
    return torch.mean(weighted_var, dim=-1, keepdim=True)

max_stars_repo_name: tqxli/dannce-pytorch

max_issues_repo_name: tqxli/dannce-pytorch

max_forks_repo_name: tqxli/dannce-pytorch


text:
#!/usr/bin/env python
# -*- coding:utf-8 -*-
import warnings
from typing import Dict, Any, Iterable, Optional, Set, Union, Tuple, List

import numpy as np
import torch
from nltk.corpus import wordnet as wn

from .wsd_baseline import MostFrequentSenseWSDTaskEvaluator, numeric
from dataset import WSDTaskDataset
from dataset.gloss_embeddings import SREFLemmaEmbeddingsDataset
from dataset.utils import tensor_to_numpy, numpy_to_tensor, batch_tile
from model.similarity import CosineSimilarity, DotProductSimilarity

class FrozenBERTKNNWSDTaskEvaluator(MostFrequentSenseWSDTaskEvaluator):

    def __init__(self,
                 evaluation_dataset: WSDTaskDataset,
                 lemma_key_embeddings_dataset: SREFLemmaEmbeddingsDataset,
                 gloss_projection_head: torch.nn.Module,
                 context_projection_head: Optional[torch.nn.Module] = None,
                 target_pos: Tuple[str] = ("n","v","a","s","r"),
                 similarity_module: Union[str, torch.nn.Module] = "cosine",
                 evaluation_category: str = "lemma",
                 ground_truth_lemma_keys_field_name: str = "ground_truth_lemma_keys",
                 entity_embedding_field_name: str = "entity_span_avg_vectors",
                 breakdown_attributes: Optional[Iterable[Set[str]]] = None,
                 device: Optional[str] = "cpu",
                 verbose: bool = False,
                 **kwargs_dataloader):
        super().__init__(
            evaluation_dataset=evaluation_dataset,
            evaluation_category=evaluation_category,
            ground_truth_lemma_keys_field_name=ground_truth_lemma_keys_field_name,
            breakdown_attributes=breakdown_attributes,
            device=device,
            verbose=verbose,
            **kwargs_dataloader)

        self._entity_embedding_field_name = entity_embedding_field_name
        self._lemma_key_embeddings_dataset = lemma_key_embeddings_dataset
        self._target_pos = target_pos
        self._gloss_prjection_head = gloss_projection_head
        if context_projection_head is None:
            warnings.warn(f"gloss_projection_head will be used as context_projection_head.")
            self._context_projection_head = gloss_projection_head
        else:
            self._context_projection_head = context_projection_head

        if isinstance(similarity_module, str):
            if similarity_module == "cosine":
                self._similarity_module = CosineSimilarity()
            elif similarity_module == "dot":
                self._similarity_module = DotProductSimilarity()
            else:
                raise ValueError(f"unknown `similarity_module` name: {similarity_module}")
        else:
            self._similarity_module = similarity_module

    def get_lemma_key_embedding(self, lemma_key: str) -> np.ndarray:
        """
        get lemma key embedding which is precomputed using WordNet gloss corpus.

        Args:
            lemma_key: lemma key.

        Returns: lemma key embedding. shape: (n_dim, )

        """
        lst_records = self._lemma_key_embeddings_dataset.get_records_by_lemma_key(lemma_key=lemma_key)
        assert len(lst_records) == 1, f"record must be unique for each lemma key: {lemma_key}"
        record = lst_records[0]
        assert lemma_key in record["ground_truth_lemma_keys"], f"lemma key lookup failure: {lemma_key}"

        return record["embeddings"]

    def get_lemma_key_embeddings(self, lst_lemma_keys: List[str]) -> torch.Tensor:
        it = map(self.get_lemma_key_embedding, lst_lemma_keys)
        v_embs = np.stack(list(it))
        return numpy_to_tensor(v_embs)

    def predict(self, input: Dict[str, Any],
                mfs_reorder_by_lemma_counts: bool = False,
                output_tie_lemma: bool = False):
        """
        predict the most plausible sense based on conditional probability

        @param input:
        @param use_generated_code_probability:
        @param apply_one_hot_encoding: convert from continuous relaxed repr. to one-hot repr.
        @return:
        """

        # when ties happen, then we fall back to most frequent sense.
        ties_fallback_to_mfs = (output_tie_lemma == False)

        lemma_name = input["lemma"]
        pos = input["pos"]

        # if unsupported part-of-speech tag, then fall back to most frequent sense method.
        if pos not in self._target_pos:
            return super().predict(input, reorder_by_lemma_counts=mfs_reorder_by_lemma_counts, output_tie_lemma=output_tie_lemma)

        # get candidates
        lst_candidate_lemmas = self.get_candidate_lemmas_from_wordnet(lemma_name, pos)
        n_candidates = len(lst_candidate_lemmas)
        lst_candidate_lemma_keys = [lemma.key() for lemma in lst_candidate_lemmas]

        # lookup candidate lemma key embeddings from WordNet gloss embeddings dataset
        device = input["entity_embeddings"].device
        t_candidate_lemma_key_embeddings = self.get_lemma_key_embeddings(lst_candidate_lemma_keys).to(device)

        # query embedding:
        # shape: (1, n_dim) if entity_embedding_field_name = "entity_span_avg_vectors"
        # shape: (1, n_subwords, n_dim) if entity_embedding_field_name = "entity_span_avg_vectors"
        t_query_embedding = input[self._entity_embedding_field_name]
        if t_query_embedding.ndim == 3:
            # average along subword dimension (=2nd dimension)
            t_query_embedding = torch.mean(t_query_embedding, dim=1)

        # inference using k-NN method.
        t_candidate_lemma_key_embeddings = self._gloss_prjection_head.predict(t_candidate_lemma_key_embeddings, is_gloss_embeddings=True)
        t_query_embedding = self._context_projection_head.predict(t_query_embedding, is_gloss_embeddings=False)
        t_query_embeddings = torch.tile(t_query_embedding, (n_candidates, 1))
        assert t_candidate_lemma_key_embeddings.shape == t_query_embeddings.shape, f"query and candidate shape mismatch."

        # calculate similarity
        t_sim_score = self._similarity_module(t_query_embeddings, t_candidate_lemma_key_embeddings)
        lst_metric_scores = tensor_to_numpy(t_sim_score).tolist()

        # return top-k lemma keys
        if ties_fallback_to_mfs:
            lst_sense_freq_scores = self.score_by_sense_frequency(lst_lemmas=lst_candidate_lemmas, reorder_by_lemma_count=mfs_reorder_by_lemma_counts)
            # order by metric first, then order by wordnet frequency rank.
            lst_scores = (lst_metric_scores, lst_sense_freq_scores)
        else:
            lst_scores = lst_metric_scores

        lst_predictions = self.return_top_k_lemma_keys(lst_candidate_lemmas, lst_scores, multiple_output=output_tie_lemma)
        dict_prediction_scores = dict(zip(lst_candidate_lemma_keys, lst_metric_scores))
        return lst_predictions, dict_prediction_scores

    def __iter__(self):
        it = super().__iter__()
        for inputs_for_predictor, inputs_for_evaluator, ground_truthes, predictions, prediction_scores, dict_metrics in it:
            if self.verbose:
                print(f"ground truth: {inputs_for_evaluator['lemma']}|{inputs_for_predictor['pos']}")
                lexnames = list(map(self._lemma_to_lexname, ground_truthes))
                ground_truth_synset_ids = list(map(self._lemma_to_synset_id, ground_truthes))
                for lexname, synset_id, ground_truth in zip(lexnames, ground_truth_synset_ids, ground_truthes):
                    print(f"\t{lexname}-{synset_id}-{ground_truth}")
                print("-------------")
            yield inputs_for_predictor, inputs_for_evaluator, ground_truthes, predictions, prediction_scores, dict_metrics

    def _print_verbose(self, lst_tup_lemma_and_scores: List[Tuple[wn.lemma, Union[numeric, Tuple[numeric]]]]):
        print(f"similarity_module: {self._similarity_module.__class__.__name__}")
        print(f"candidates:")
        for lemma, scores in lst_tup_lemma_and_scores:
            lexname = lemma.synset().lexname()
            synset_id = lemma.synset().name()
            if isinstance(scores, float):
                print(f"\t{lexname}-{synset_id}-{lemma.key()}: {scores:1.6f}")
            else:
                print(f"\t{lexname}-{synset_id}-{lemma.key()}: {scores}")


text:
# ===================================================
#  UMAP Fit and Transform Operations Test cases
#  (not really fitting anywhere else)
# ===================================================

from sklearn.datasets import make_blobs
from sklearn.cluster import KMeans
from sklearn.metrics import adjusted_rand_score, pairwise_distances
from sklearn.preprocessing import normalize
from nose.tools import assert_equal, assert_less, assert_raises
from numpy.testing import assert_array_equal
from umap import UMAP
from umap.spectral import component_layout
import numpy as np
import scipy.sparse
import pytest
import warnings
from umap.distances import pairwise_special_metric
from umap.utils import disconnected_vertices
from scipy.sparse import csr_matrix

# Transform isn't stable under batching; hard to opt out of this.
# @SkipTest
# def test_scikit_learn_compatibility():
#     check_estimator(UMAP)

# This test is currently to expensive to run when turning
# off numba JITting to detect coverage.
# @SkipTest
# def test_umap_regression_supervision(): # pragma: no cover
#     boston = load_boston()
#     data = boston.data
#     embedding = UMAP(n_neighbors=10,
#                      min_dist=0.01,
#                      target_metric='euclidean',
#                      random_state=42).fit_transform(data, boston.target)
#

# Umap Clusterability
def test_blobs_cluster():
    data, labels = make_blobs(n_samples=500, n_features=10, centers=5)
    embedding = UMAP(n_epochs=100).fit_transform(data)
    assert_equal(adjusted_rand_score(labels, KMeans(5).fit_predict(embedding)), 1.0)

# Multi-components Layout
def test_multi_component_layout():
    data, labels = make_blobs(
        100, 2, centers=5, cluster_std=0.5, center_box=[-20, 20], random_state=42
    )

    true_centroids = np.empty((labels.max() + 1, data.shape[1]), dtype=np.float64)

    for label in range(labels.max() + 1):
        true_centroids[label] = data[labels == label].mean(axis=0)

    true_centroids = normalize(true_centroids, norm="l2")

    embedding = UMAP(n_neighbors=4, n_epochs=100).fit_transform(data)
    embed_centroids = np.empty((labels.max() + 1, data.shape[1]), dtype=np.float64)
    embed_labels = KMeans(n_clusters=5).fit_predict(embedding)

    for label in range(embed_labels.max() + 1):
        embed_centroids[label] = data[embed_labels == label].mean(axis=0)

    embed_centroids = normalize(embed_centroids, norm="l2")

    error = np.sum((true_centroids - embed_centroids) ** 2)

    assert_less(error, 15.0, msg="Multi component embedding to far astray")

# Multi-components Layout
def test_multi_component_layout_precomputed():
    data, labels = make_blobs(
        100, 2, centers=5, cluster_std=0.5, center_box=[-20, 20], random_state=42
    )
    dmat = pairwise_distances(data)

    true_centroids = np.empty((labels.max() + 1, data.shape[1]), dtype=np.float64)

    for label in range(labels.max() + 1):
        true_centroids[label] = data[labels == label].mean(axis=0)

    true_centroids = normalize(true_centroids, norm="l2")

    embedding = UMAP(n_neighbors=4, metric="precomputed", n_epochs=100).fit_transform(
        dmat
    )
    embed_centroids = np.empty((labels.max() + 1, data.shape[1]), dtype=np.float64)
    embed_labels = KMeans(n_clusters=5).fit_predict(embedding)

    for label in range(embed_labels.max() + 1):
        embed_centroids[label] = data[embed_labels == label].mean(axis=0)

    embed_centroids = normalize(embed_centroids, norm="l2")

    error = np.sum((true_centroids - embed_centroids) ** 2)

    assert_less(error, 15.0, msg="Multi component embedding to far astray")

@pytest.mark.parametrize("num_isolates", [1, 5])
@pytest.mark.parametrize("metric", ["jaccard", "hellinger", "cosine"])
@pytest.mark.parametrize("force_approximation", [True, False])
def test_disconnected_data(num_isolates, metric, force_approximation):
    disconnected_data = np.random.choice(
        a=[False, True], size=(10, 20), p=[0.66, 1 - 0.66]
    )
    # Add some disconnected data for the corner case test
    disconnected_data = np.vstack(
        [disconnected_data, np.zeros((num_isolates, 20), dtype="bool")]
    )
    new_columns = np.zeros((num_isolates + 10, num_isolates), dtype="bool")
    for i in range(num_isolates):
        new_columns[10 + i, i] = True
    disconnected_data = np.hstack([disconnected_data, new_columns])

    with warnings.catch_warnings(record=True) as w:
        model = UMAP(
            n_neighbors=3,
            metric=metric,
            force_approximation_algorithm=force_approximation,
        ).fit(disconnected_data)
        assert len(w) >= 1  # at least one warning should be raised here
        # we can't guarantee the order that the warnings will be raised in so check them all.
        flag = 0
        if num_isolates == 1:
            warning_contains = "A few of your vertices"
        elif num_isolates > 1:
            warning_contains = "A large number of your vertices"
        for i in range(len(w)):
            flag += warning_contains in str(w[i].message)
        assert flag == 1
        # Check that the first isolate has no edges in our umap.graph_
        isolated_vertices = disconnected_vertices(model)
        assert isolated_vertices[10] == True
        number_of_nan = np.sum(np.isnan(model.embedding_[isolated_vertices]))
        assert number_of_nan >= num_isolates*model.n_components

@pytest.mark.parametrize("num_isolates", [1])
@pytest.mark.parametrize("sparse", [True, False])
def test_disconnected_data_precomputed(num_isolates, sparse):
    disconnected_data = np.random.choice(
        a=[False, True], size=(10, 20), p=[0.66, 1 - 0.66]
    )
    # Add some disconnected data for the corner case test
    disconnected_data = np.vstack(
        [disconnected_data, np.zeros((num_isolates, 20), dtype="bool")]
    )
    new_columns = np.zeros((num_isolates + 10, num_isolates), dtype="bool")
    for i in range(num_isolates):
        new_columns[10 + i, i] = True
    disconnected_data = np.hstack([disconnected_data, new_columns])
    dmat = pairwise_special_metric(disconnected_data)
    if sparse:
        dmat = csr_matrix(dmat)
    model = UMAP(n_neighbors=3, metric="precomputed", disconnection_distance=1).fit(
        dmat
    )

    # Check that the first isolate has no edges in our umap.graph_
    isolated_vertices = disconnected_vertices(model)
    assert isolated_vertices[10] == True
    number_of_nan = np.sum(np.isnan(model.embedding_[isolated_vertices]))
    assert number_of_nan >= num_isolates * model.n_components

# ---------------
# Umap Transform
# --------------

def test_bad_transform_data(nn_data):
    u = UMAP().fit([[1, 1, 1, 1]])
    assert_raises(ValueError, u.transform, [[0, 0, 0, 0]])

# Transform Stability
# -------------------
def test_umap_transform_embedding_stability(iris, iris_subset_model, iris_selection):
    """Test that transforming data does not alter the learned embeddings

    Issue #217 describes how using transform to embed new data using a
    trained UMAP transformer causes the fitting embedding matrix to change
    in cases when the new data has the same number of rows as the original
    training data.
    """

    data = iris.data[iris_selection]
    fitter = iris_subset_model
    original_embedding = fitter.embedding_.copy()

    # The important point is that the new data has the same number of rows
    # as the original fit data
    new_data = np.random.random(data.shape)
    _ = fitter.transform(new_data)

    assert_array_equal(
        original_embedding,
        fitter.embedding_,
        "Transforming new data changed the original embeddings",
    )

    # Example from issue #217
    a = np.random.random((100, 10))
    b = np.random.random((100, 5))

    umap = UMAP(n_epochs=100)
    u1 = umap.fit_transform(a[:, :5])
    u1_orig = u1.copy()
    assert_array_equal(u1_orig, umap.embedding_)

    _ = umap.transform(b)
    assert_array_equal(u1_orig, umap.embedding_)

# -----------
# UMAP Update
# -----------
def test_umap_update(iris, iris_subset_model, iris_selection, iris_model):

    new_data = iris.data[~iris_selection]
    new_model = iris_subset_model
    new_model.update(new_data)

    comparison_graph = scipy.sparse.vstack(
        [iris_model.graph_[iris_selection], iris_model.graph_[~iris_selection]]
    )
    comparison_graph = scipy.sparse.hstack(
        [comparison_graph[:, iris_selection], comparison_graph[:, ~iris_selection]]
    )

    error = np.sum(np.abs((new_model.graph_ - comparison_graph).data))

    assert_less(error, 1.0)

# -----------------
# UMAP Graph output
# -----------------
def test_umap_graph_layout():
    data, labels = make_blobs(n_samples=500, n_features=10, centers=5)
    model = UMAP(n_epochs=100, transform_mode="graph")
    graph = model.fit_transform(data)
    assert scipy.sparse.issparse(graph)
    nc, cl = scipy.sparse.csgraph.connected_components(graph)
    assert_equal(nc, 5)

    new_graph = model.transform(data[:10] + np.random.normal(0.0, 0.1, size=(10, 10)))
    assert scipy.sparse.issparse(graph)
    assert new_graph.shape[0] == 10

# ------------------------
# Component layout options
# ------------------------

def test_component_layout_options(nn_data):
    dmat = pairwise_distances(nn_data[:1000])
    n_components = 5
    component_labels = np.repeat(np.arange(5), dmat.shape[0] // 5)
    single = component_layout(
        dmat,
        n_components,
        component_labels,
        2,
        np.random,
        metric="precomputed",
        metric_kwds={"linkage": "single"},
    )
    average = component_layout(
        dmat,
        n_components,
        component_labels,
        2,
        np.random,
        metric="precomputed",
        metric_kwds={"linkage": "average"},
    )
    complete = component_layout(
        dmat,
        n_components,
        component_labels,
        2,
        np.random,
        metric="precomputed",
        metric_kwds={"linkage": "complete"},
    )

    assert single.shape[0] == 5
    assert average.shape[0] == 5
    assert complete.shape[0] == 5

    assert not np.all(single == average)
    assert not np.all(single == complete)
    assert not np.all(average == complete)

max_stars_repo_name: nidhalrahali/umap

max_issues_repo_name: nidhalrahali/umap

max_forks_repo_name: nidhalrahali/umap

reason: import numpy,from numpy,import scipy,from scipy


text:
import numpy as np
import math

from model.alphabet import Alphabet

def divisor_list (m):
    divisor_ls = []
    for i in range(2, m):
        if m % i == 0:
            divisor_ls.append(i)

    return divisor_ls

# alphabet = Alphabet("ABCDEFGHIJKLMNOPQRSTUVWXYZ")
# matrix = np.array([[5, 6], [2, 3]])

def is_invertible(matrix, alphabet):

    det = round(np.linalg.det(matrix))
    if det == 0: return False
    divisor_ls = divisor_list(alphabet.alpha_size)
    for d in divisor_ls:
        if det % d == 0: return False

    return True
#
# print(is_invertible(matrix, alphabet))
def has_valid_length (n):
    if (n == 1):
        return True, 1

    for i in range(2, math.ceil(n/2 +1)):
        if i*i == n:
            return True, i

    return False, -1

# check invertible matrix mod 26

# tao matrix tu password:
# for numerical password
def to_matrix(password, n, type_password, alphabet):
    if type_password == 0:
        matrix = np.array([int(digit) for digit in password]).reshape((n, n))

    # for letter password:
    else:
        matrix = np.array([alphabet.let_map_num.get(letter) for letter in password]).reshape((n, n))

    return matrix


text:
import numpy as np
import os
from os.path import sep, dirname, realpath
import pyqtgraph as pg
from pyqtgraph import opengl as gl
from pyqtgraph.Qt import QtCore, QtGui
import sys
import scipy.optimize as optimize
import matplotlib.pyplot as plt

path = os.path.dirname(os.path.abspath(__file__))
path = path.split("throttle")[0]
sys.path.insert(0, path)
import analisis_functions as af  # They are good AF

class Visualizer(object):
    def __init__(self):

        find_braking_values = False
        find_coast_values = False
        find_acceleration_values = False
        # Get the directory for the data
        data_dir = dirname(realpath(__file__)) + r"\data\regular"
        data_file_name = [
            "full_data_path.npy",
            "#acc_data_path.npy",
            "#break_data_path.npy",
            "#coast_data_path.npy",
        ]
        # Initiate the data array
        full_data = np.empty((3,))
        # iterating through all files and adding them to a huge array
        for filename in os.listdir(data_dir):
            if filename not in data_file_name:
                continue
            filepath = data_dir + sep + filename
            # Loading data
            data = np.load(filepath)
            full_data = np.vstack((full_data, data))
        self.full_data = full_data

        # get the acceleration
        self.full_data_grad = np.gradient(self.full_data[:, 2], self.full_data[:, 1])
        self.full_data_grad = np.nan_to_num(
            self.full_data_grad, nan=0, posinf=0, neginf=0
        )
        self.full_data_grad = np.hstack(
            (
                np.transpose(self.full_data[:, 0][np.newaxis]),
                np.transpose(self.full_data[:, 2][np.newaxis]),
                np.transpose(self.full_data_grad[np.newaxis]),
            )
        )

        if find_braking_values:
            acceleration = abs(self.full_data_grad[:, 2])
            acceleration = np.sort(acceleration)
            acceleration = np.around(acceleration, decimals=3)
            normal_mask = (acceleration > 3498) & (acceleration < 3502)
            boost_mask = (
                (self.full_data_grad[:, 0] > 1)
                & (acceleration > 4491)
                & (acceleration < 4492)
            )
            normal_acceleration = acceleration[normal_mask]
            boost_acceleration = acceleration[boost_mask]
            print(f"Normal braking force = {np.median(normal_acceleration)}")
            print(f"Boost braking force = {np.median(boost_acceleration)}")
            print(
                f"Boost braking dif = {np.median(boost_acceleration)-np.median(normal_acceleration)}"
            )

        if find_coast_values:
            acceleration = abs(self.full_data_grad[:, 2])
            acceleration = np.sort(acceleration)
            acceleration = np.around(acceleration, decimals=4)
            normal_mask = (acceleration > 527.661) & (acceleration < 527.663)
            normal_acceleration = acceleration[normal_mask]
            hist = np.histogram(normal_acceleration, bins=100)
            print(f"Coasting force = {np.median(normal_acceleration)}")

        if find_acceleration_values:
            acceleration = abs(self.full_data_grad[:, 2])
            acceleration = np.around(acceleration, decimals=5)
            normal_mask = (
                (acceleration > 0)
                & (acceleration < 3500)
                & (self.full_data_grad[:, 0] <= 1)
            )
            normal_acceleration = acceleration[normal_mask]
            normal_acceleration = np.sort(normal_acceleration)

            stable_acceleration = acceleration / self.full_data_grad[:, 0]
            stable_func_grad = np.gradient(
                stable_acceleration, self.full_data_grad[:, 1]
            )
            stable_func_grad = np.nan_to_num(
                stable_func_grad, nan=0, posinf=0, neginf=0
            )

            close_to_transition_mask = (
                (abs(self.full_data_grad[:, 1]) > 1390)
                & (abs(self.full_data_grad[:, 1]) < 1500)
                & (abs(stable_func_grad) > 50)
            )
            speed_transition = self.full_data_grad[close_to_transition_mask, 1]

            plt.scatter(
                abs(self.full_data_grad[:, 1]), abs(stable_acceleration), marker="."
            )
            plt.title("Scatter plot pythonspot.com")
            plt.xlabel("x")
            plt.ylabel("y")
            # plt.show()

            print(
                f"Max acceleration force = {(np.max(acceleration)-990.413) / 0.99609375}"
            )
            print(f"Normal acceleration force = {np.median(normal_acceleration)}")
            print(f"Acceleration transition (speed) = {np.min(abs(speed_transition))}")
            print(f"Max normal speed = {np.max(abs(speed_transition))}")
            print(f"Max boost speed = {np.max(abs(self.full_data_grad[:, 1]))}")

        def array_to_onebyone(array1, array2, function):
            if array1.size != array2.size:
                raise ValueError(
                    "all the input arrays must have same number of dimensions"
                )
            total_result = np.empty((0,))
            for index in range(array1.size):
                result = function(array1[index], array2[index])
                total_result = np.hstack((total_result, result))
            return total_result

        def tick_acceleration(
            throttle,
            speed,
            BRAKING_ACCELERATION=3501.255,
            BOOST_IMPULSE=990.413,
            COASTING_FORCE=527.6618,
            MAX_NORMAL_ACCELERATION=1610.733,
            ACCELERATION_TRANSITION_POINT=[1400.55, 160],
            MAX_NORMAL_SPEED=1409.99,
            MAX_BOOST_SPEED=2300.001,
        ):

            coast_limit = 0.0117647058823529
            sign_throttle = np.sign(throttle)
            sign_speed = np.sign(speed)
            fase = sign_speed * sign_throttle
            if abs(throttle) < coast_limit:
                mode = "coast"
            elif fase > 0:
                mode = "accelerating"
            else:
                mode = "braking"

            if abs(throttle) > 1:
                boost = True
            else:
                boost = False
            if mode == "braking":
                current_acceleration = BRAKING_ACCELERATION * sign_throttle
            if mode == "coast":
                current_acceleration = COASTING_FORCE * -1 * sign_speed
            if mode == "accelerating":
                if abs(speed) < ACCELERATION_TRANSITION_POINT[0]:
                    gradient = (
                        ACCELERATION_TRANSITION_POINT[1] - MAX_NORMAL_ACCELERATION
                    ) / (ACCELERATION_TRANSITION_POINT[0])

                    current_acceleration = (
                        (gradient * speed + MAX_NORMAL_ACCELERATION * sign_speed)
                        * throttle
                        * sign_speed
                    )

                elif abs(speed) < MAX_NORMAL_SPEED:
                    gradient = (-ACCELERATION_TRANSITION_POINT[1]) / (
                        MAX_NORMAL_SPEED - ACCELERATION_TRANSITION_POINT[0]
                    )

                    current_acceleration = (
                        (gradient * speed + 23898.1355932202 * sign_speed)
                        * throttle
                        * sign_speed
                    )
                else:
                    gradient = 0
                    current_acceleration = 0

            if boost:
                if speed >= MAX_BOOST_SPEED:
                    current_acceleration = 0
                else:
                    current_acceleration = tick_acceleration(
                        1, speed
                    ) + BOOST_IMPULSE * np.where(
                        np.sign(current_acceleration) == 0,
                        1,
                        np.sign(current_acceleration),
                    )
            return current_acceleration

        self.full_data_grad_func = np.hstack(
            (
                np.transpose(self.full_data_grad[:, 0][np.newaxis]),
                np.transpose(self.full_data_grad[:, 1][np.newaxis]),
                np.transpose(
                    array_to_onebyone(
                        self.full_data_grad[:, 0],
                        self.full_data_grad[:, 1],
                        tick_acceleration,
                    )[np.newaxis]
                ),
            )
        )

        # PyQtGraph stuff
        self.app = QtGui.QApplication(sys.argv)
        self.w = gl.GLViewWidget()
        self.w.opts["distance"] = 1000
        self.w.setWindowTitle("RL Throttle Visualized")
        self.w.setGeometry(0, 30, 1920, 1080)
        self.w.show()

        # Adding a grid
        gz = gl.GLGridItem(QtGui.QVector3D(1000, 1500, 1))
        gz.translate(500, 1500 / 2, 0)
        # self.w.addItem(gz)

        # Create the scatter plot
        blue_color = [0.1, 0.5, 1, 0.65]
        red_color = [1, 0.3, 0.05, 0.65]
        white_color = [1.0, 1.0, 1.0, 0.6]

        self.full_data_grad[:, 0] = self.full_data_grad[:, 0] * 1000
        scatter = gl.GLScatterPlotItem(
            pos=self.full_data_grad, size=1, color=blue_color
        )
        self.w.addItem(scatter)

        self.full_data_grad_func[:, 0] = self.full_data_grad_func[:, 0] * 1000
        scatter = gl.GLScatterPlotItem(
            pos=self.full_data_grad_func, size=1, color=red_color
        )
        self.w.addItem(scatter)

        # https://stackoverflow.com/questions/56890547/how-to-add-axis-features-labels-ticks-values-to-a-3d-plot-with-glviewwidget
        axis = af.Custom3DAxis(self.w, color=[1.0, 1.0, 1.0, 1.0])
        axis.setSize(x=1100, y=1600, z=3800)
        # Add axes labels
        axis.add_labels()
        # Add axes tick values
        axis.add_tick_values(
            xticks=list(range(0, 11, 1)) + ["Throttle"],
            yticks=list(range(0, 1600, 100)) + ["Current Speed"],
            zticks=list(range(0, 3800, 100)) + ["Acceleration"],
        )
        self.w.addItem(axis)

    def start(self):
        if (sys.flags.interactive != 1) or not hasattr(QtCore, "PYQT_VERSION"):
            QtGui.QApplication.instance().exec_()

    def update(self):
        pass

    def animation(self):
        timer = QtCore.QTimer()
        timer.timeout.connect(self.update)
        timer.start(20)
        self.start()

if __name__ == "__main__":
    v = Visualizer()
    v.animation()

max_stars_repo_name: skyborgff/Unraveling-Rocket-League

max_issues_repo_name: skyborgff/Unraveling-Rocket-League

max_forks_repo_name: skyborgff/Unraveling-Rocket-League

reason: import numpy,import scipy


text:
from latexqol.conversion import copy_table
import numpy as np

head = ["head1", "head2", "head3", "head4"]
col1 = np.arange(0, 20, 1)
col2 = np.arange(20, 40, 1)
col3 = np.arange(40, 60, 1)
col4 = np.arange(60, 80, 1)

copy_table(head, col1, col2, col3, col4)

max_stars_repo_name: konradgj/latexqol

max_issues_repo_name: konradgj/latexqol

max_forks_repo_name: konradgj/latexqol


text:
import pandas as pd
import numpy as np
import os
from sklearn.model_selection import train_test_split
from model_tensorflow import train, predict

frame = "tensorflow"

class Config:
	# feature_columns = list(range(0,8))
	# label_columns = [5,6,7]
	feature_columns = list([2,5])#comment yqy
	# feature_columns = list([2]) #add yqy
	label_columns = [5]
	feature_and_label_columns = feature_columns + label_columns
	label_in_feature_columns = (lambda x, y: [x.index(i) for i in y])(feature_columns, label_columns)

	predict_day = 1

	# input_size = len(feature_columns)#comment yqy
	input_size = len( list([2]))#add yqy
	output_size = len(label_columns)

	hidden_size = 128
	lstm_layers = 2
	dropout_rate = 0.2
	time_step = 5

	# do_train = True
	do_train = True
	do_predict = True
	add_train = False
	shuffle_train_data = True

	# train_data_rate = 0.95 #comment yqy
	train_data_rate = 1 #add yqy
	valid_data_rate = 0.15

	batch_size = 64
	learning_rate = 0.001
	epoch = 20
	patience = 5
	random_seed = 42

	do_continue_train = False
	continue_flag = ""
	if do_continue_train:
		shuffle_train_data = False
		batch_size = 1
		continue_flag = "continue_"

	#comment yqy
	train_data_path = "./data/stock_data.csv"
	model_save_path = "./checkpoint/"
	figure_save_path = "./figure/"
	#comment end
	# add yqy
	# train_data_path = "./data/stock_data_30.csv"
	# model_save_path = "./checkpoint/30/"
	# figure_save_path = "./figure/30/"
	# add end
	do_figure_save = False
	if not os.path.exists(model_save_path):
		os.mkdir(model_save_path)
	if not os.path.exists(figure_save_path):
		os.mkdir(figure_save_path)

	used_frame = frame
	model_postfix = {"pytorch": ".pth", "keras": ".h5", "tensorflow": ".ckpt"}
	model_name = "model_" + continue_flag + used_frame + model_postfix[used_frame]

class Data:
	def __init__(self, config):
		self.config = config
		self.data, self.data_column_name = self.read_data()

		self.data_num = self.data.shape[0]
		self.train_num = int(self.data_num * self.config.train_data_rate)

		self.mean = np.mean(self.data, axis=0)
		self.std = np.std(self.data, axis=0)
		self.norm_data = (self.data - self.mean) / self.std

		self.start_num_in_test = 0

	def read_data(self):
		init_data = pd.read_csv(self.config.train_data_path,
		                        usecols=self.config.feature_and_label_columns)
		return init_data.values, init_data.columns.tolist()

	def get_train_and_valid_data(self):
		# feature_data = self.norm_data[:self.train_num] # comment yqy
		feature_data = self.norm_data[:self.train_num][:,1][:,np.newaxis] # add yqy
		label_data = self.norm_data[self.config.predict_day: self.config.predict_day + self.train_num,
		             self.config.label_in_feature_columns]
		if not self.config.do_continue_train:
			train_x = [feature_data[i:i + self.config.time_step] for i in range(self.train_num - self.config.time_step)]
			train_y = [label_data[i:i + self.config.time_step] for i in range(self.train_num - self.config.time_step)]
		else:
			train_x = [
				feature_data[start_index + i * self.config.time_step: start_index + (i + 1) * self.config.time_step]
				for start_index in range(self.config.time_step)
				for i in range((self.train_num - start_index) // self.config.time_step)]
			train_y = [
				label_data[start_index + i * self.config.time_step: start_index + (i + 1) * self.config.time_step]
				for start_index in range(self.config.time_step)
				for i in range((self.train_num - start_index) // self.config.time_step)]

		train_x, train_y = np.array(train_x), np.array(train_y)

		train_x, valid_x, train_y, valid_y = train_test_split(train_x, train_y, test_size=self.config.valid_data_rate,
		                                                      random_state=self.config.random_seed,
		                                                      shuffle=self.config.shuffle_train_data)
		return train_x, valid_x, train_y, valid_y

	def get_test_data(self, return_label_data=False):
		feature_data = self.norm_data[self.train_num:]
		self.start_num_in_test = feature_data.shape[0] % self.config.time_step
		time_step_size = feature_data.shape[0] // self.config.time_step

		test_x = [feature_data[self.start_num_in_test + i * self.config.time_step: self.start_num_in_test + (
				i + 1) * self.config.time_step]
		          for i in range(time_step_size)]
		if return_label_data:
			label_data = self.norm_data[self.train_num + self.start_num_in_test:, self.config.label_in_feature_columns]
			return np.array(test_x), label_data
		return np.array(test_x)

	# add yqy
	def get_test_data_yqy(self, test_data_yqy=None):
		if test_data_yqy is None:
			test_data_yqy = []
		# test_data_yqy=test_data_yqy[1:21]
		feature_data=(test_data_yqy - self.mean) / self.std
		test_x=[feature_data]
		return np.array(test_x)
	# add end

def draw(config, origin_data, predict_norm_data):
	label_norm_data = origin_data.norm_data[origin_data.train_num + origin_data.start_num_in_test:,
	                  config.label_in_feature_columns]
	assert label_norm_data.shape[0] == predict_norm_data.shape[
		0], "The element number in origin and predicted data is different"

	label_name = [origin_data.data_column_name[i] for i in config.label_in_feature_columns]
	label_column_num = len(config.label_columns)

	loss = np.mean((label_norm_data[config.predict_day:] - predict_norm_data[:-config.predict_day]) ** 2, axis=0)
	print("The mean squared error of stock {} is ".format(label_name), loss)

	label_X = range(origin_data.data_num - origin_data.train_num - origin_data.start_num_in_test)
	predict_X = [x + config.predict_day for x in label_X]

	label_data = label_norm_data * origin_data.std[config.label_in_feature_columns] + \
	             origin_data.mean[config.label_in_feature_columns]

	predict_data = predict_norm_data * origin_data.std[config.label_in_feature_columns] + \
	               origin_data.mean[config.label_in_feature_columns]

	print(label_data)
	print("____________________________________________")
	print(predict_data)

def draw_yqy(config, origin_data, predict_norm_data,mean_yqy,std_yqy):# 这里origin_data等同于test_data_values_yqy
	label_norm_data = (origin_data - mean_yqy) / std_yqy
	assert label_norm_data.shape[0] == predict_norm_data.shape[0], "The element number in origin and predicted data is different"

	#label_norm_data=label_norm_data[:,1]
	label_name = 'high'
	label_column_num = 1

	loss = np.mean((label_norm_data[config.predict_day:][:,1][:,np.newaxis] - predict_norm_data[:-config.predict_day][0:]) ** 2, axis=0)

	# loss = np.mean((label_norm_data[config.predict_day:][:,5][:,np.newaxis] - predict_norm_data[:-config.predict_day][0:]) ** 2, axis=0)
	# loss2 = np.mean((label_norm_data[config.predict_day:][:,6][:,np.newaxis] - predict_norm_data[:-config.predict_day][0:]) ** 2, axis=0)
	# loss3 = np.mean((label_norm_data[config.predict_day:][:,7][:,np.newaxis] - predict_norm_data[:-config.predict_day][0:]) ** 2, axis=0)

	print("The mean squared error of stock {} is ".format(label_name), loss)

	# label_X = range(origin_data.data_num - origin_data.train_num - origin_data.start_num_in_test)
	# predict_X = [x + config.predict_day for x in label_X]

	label_data = label_norm_data[:,1] * std_yqy[1]+ mean_yqy[1]

	predict_data = predict_norm_data * std_yqy[1]+ mean_yqy[1]

	print(label_data)
	print(predict_data)
	# print(label_data[-1])
	# print(predict_data[-1][0])

def main(config):
	np.random.seed(config.random_seed)
	data_gainer = Data(config)

	# add yqy
	mean_yqy=Data(config).mean
	std_yqy=Data(config).std
	#add end

	if config.do_train:
		train_X, valid_X, train_Y, valid_Y = data_gainer.get_train_and_valid_data()
		train(config, train_X, train_Y, valid_X, valid_Y)

	if config.do_predict:
		# add yqy
		test_data_yqy = pd.read_csv("./data/test_data.csv",usecols=list([2,5]))
		test_data_values_yqy=test_data_yqy.values
		# test_data_yqy=[104.3,104.39]
		test_X =data_gainer.get_test_data_yqy(test_data_values_yqy)
		# add end
		# test_X, test_Y = data_gainer.get_test_data(return_label_data=True)# comment yqy
		# pred_result = predict(config, test_X)
		pred_result = predict(config,test_X[:,:,0][:,:,np.newaxis])
		# draw(config, data_gainer, pred_result)# comment yqy
		draw_yqy(config, test_data_values_yqy, pred_result,mean_yqy,std_yqy)

if __name__ == "__main__":
	import argparse

	parser = argparse.ArgumentParser()
	# parser.add_argument("-t", "--do_train", default=False, type=bool, help="whether to train")
	# parser.add_argument("-p", "--do_predict", default=True, type=bool, help="whether to train")
	# parser.add_argument("-b", "--batch_size", default=64, type=int, help="batch size")
	# parser.add_argument("-e", "--epoch", default=20, type=int, help="epochs num")
	args = parser.parse_args()

	con = Config()
	for key in dir(args):
		if not key.startswith("_"):
			setattr(con, key, getattr(args, key))

	main(con)

max_stars_repo_name: PasaLab/YAO

max_issues_repo_name: PasaLab/YAO

max_forks_repo_name: PasaLab/YAO


text:
"""
GAN_graph.py: Tensorflow Graph for the Generative Adversarial Networks
"""
__author__      = "Khalid M. Kahloot"
__copyright__   = "Copyright 2019, Only for professionals"

from base.base_graph import BaseGraph
import tensorflow as tf
import numpy as np
from tqdm import tqdm

from networks.dense_net import DenseNet

'''
This is the Main GANGraph.
'''

class GANGraph(BaseGraph):
    def __init__(self, network_params, sigma_act=tf.nn.softplus,
                 transfer_fct=tf.nn.relu, learning_rate=1e-4,
                 kinit=tf.contrib.layers.xavier_initializer(), batch_size=32,
                 reuse=None, dropout=0.2):
        super().__init__(learning_rate)

        self.width = network_params['input_width']
        self.height = network_params['input_height']
        self.nchannel = network_params['input_nchannels']

        self.hidden_dim = network_params['hidden_dim']
        self.latent_dim = network_params.latent_dim
        self.num_layers = network_params['num_layers']  # Num of Layers in P(x|z)
        self.l2 = network_params.l2
        self.dropout = dropout

        self.sigma_act = sigma_act  # Actfunc for NN modeling variance

        self.x_flat_dim = self.width * self.height * self.nchannel

        self.transfer_fct = transfer_fct
        self.kinit = kinit
        self.bias_init = tf.constant_initializer(0.0)
        self.batch_size = batch_size

        self.reuse = reuse

    def build_graph(self):
        self.create_inputs()
        self.create_graph()
        self.create_loss_optimizer()

    def create_inputs(self):
        with tf.variable_scope('inputs', reuse=self.reuse):
            self.x_batch = tf.placeholder(tf.float32, [self.batch_size, self.width, self.height, self.nchannel], name='x_batch')
            self.x_batch_flat = tf.reshape(self.x_batch, [-1, self.x_flat_dim])
            self.z_batch = tf.placeholder(tf.float32, [self.batch_size, self.latent_dim], name='z_batch')

    '''
    ------------------------------------------------------------------------------
                                     GRAPH FUNCTIONS
    ------------------------------------------------------------------------------
    '''

    def create_graph(self):
        print('\n[*] Defining Generator...')
        with tf.variable_scope('generator', reuse=self.reuse):
            self.Gx = self.create_decoder(input_=self.z_batch,
                                          hidden_dim=self.hidden_dim,
                                          output_dim=self.x_flat_dim,
                                          num_layers=self.num_layers,
                                          transfer_fct=self.transfer_fct,
                                          act_out=tf.nn.sigmoid,
                                          reuse=self.reuse,
                                          kinit=self.kinit,
                                          bias_init=self.bias_init,
                                          drop_rate=self.dropout,
                                          prefix='g_'
                                          )

            self.x_recons = tf.reshape(self.Gx.output, [-1, self.width, self.height, self.nchannel])
            self.x_recons_flat = tf.reshape(self.Gx.output, [-1, self.x_flat_dim])
        print('\n[*] Defining Discriminator...')
        with tf.variable_scope('discriminator', reuse=self.reuse):
            self.D_real_logit = self.create_encoder(input_=self.x_batch_flat,
                                                    hidden_dim=self.hidden_dim,
                                                    output_dim=self.latent_dim,
                                                    num_layers=self.num_layers,
                                                    transfer_fct=self.transfer_fct,
                                                    act_out=None,
                                                    reuse=self.reuse,
                                                    kinit=self.kinit,
                                                    bias_init=self.bias_init,
                                                    drop_rate=self.dropout,
                                                    prefix='d_r_')

            self.D_real_prob = tf.nn.sigmoid(self.D_real_logit.output)

            self.D_fake_logit = self.create_encoder(input_=self.x_recons_flat,
                                                    hidden_dim=self.hidden_dim,
                                                    output_dim=self.latent_dim,
                                                    num_layers=self.num_layers,
                                                    transfer_fct=self.transfer_fct,
                                                    act_out=None,
                                                    reuse=self.reuse,
                                                    kinit=self.kinit,
                                                    bias_init=self.bias_init,
                                                    drop_rate=self.dropout,
                                                    prefix='d_f_')

            self.D_fake_prob = tf.nn.sigmoid(self.D_fake_logit.output)

    '''
    ------------------------------------------------------------------------------
                                     ENCODER-DECODER
    ------------------------------------------------------------------------------
    '''

    def create_encoder(self, input_, hidden_dim, output_dim, num_layers, transfer_fct, \
                       act_out, reuse, kinit, bias_init, drop_rate, prefix):
        latent_ = DenseNet(input_=input_,
                           hidden_dim=hidden_dim,
                           output_dim=output_dim,
                           num_layers=num_layers,
                           transfer_fct=transfer_fct,
                           act_out=act_out,
                           reuse=reuse,
                           kinit=kinit,
                           bias_init=bias_init,
                           drop_rate=drop_rate,
                           prefix=prefix)
        return latent_

    def create_decoder(self, input_, hidden_dim, output_dim, num_layers, transfer_fct, \
                       act_out, reuse, kinit, bias_init, drop_rate, prefix):
        recons_ = DenseNet(input_=input_,
                           hidden_dim=hidden_dim,
                           output_dim=output_dim,
                           num_layers=num_layers,
                           transfer_fct=transfer_fct,
                           act_out=act_out,
                           reuse=reuse,
                           kinit=kinit,
                           bias_init=bias_init,
                           drop_rate=drop_rate,
                           prefix=prefix)
        return recons_

    '''
    ------------------------------------------------------------------------------
                                     LOSSES
    ------------------------------------------------------------------------------
    '''

    def create_loss_optimizer(self):
        print('[*] Defining Loss Functions and Optimizer...')
        with tf.variable_scope('d_loss', reuse=self.reuse):
            # get loss for discriminator
            self.d_loss_real = tf.reduce_mean(tf.nn.sigmoid_cross_entropy_with_logits(
                logits=self.D_real_logit.output, labels=tf.ones_like(self.D_real_prob)))

            self.d_loss_fake = tf.reduce_mean(tf.nn.sigmoid_cross_entropy_with_logits(
                logits=self.D_fake_logit.output, labels=tf.zeros_like(self.D_fake_prob)))

            self.d_loss = self.d_loss_real + self.d_loss_fake

        with tf.variable_scope('g_loss', reuse=self.reuse):
            # get loss for generator
            self.g_loss = tf.reduce_mean(
                tf.nn.sigmoid_cross_entropy_with_logits(
                    logits=self.D_fake_logit.output, labels=tf.ones_like(self.D_fake_prob)))

        with tf.name_scope('reconstruct'):
            self.reconstruction = self.get_ell(self.x_batch_flat, self.x_recons_flat)

        t_vars = tf.trainable_variables()
        d_r_vars = [var for var in t_vars if 'd_r_' in var.name]
        d_f_vars = [var for var in t_vars if 'd_f_' in var.name]
        g_vars = [var for var in t_vars if 'g_' in var.name]

        with tf.variable_scope("L2_loss", reuse=self.reuse):
            self.gL2_loss = tf.reduce_sum([tf.nn.l2_loss(v) for v in g_vars])
            self.drL2_loss = tf.reduce_sum([tf.nn.l2_loss(v) for v in d_r_vars])
            self.dfL2_loss = tf.reduce_sum([tf.nn.l2_loss(v) for v in d_f_vars])
            self.reconL2_loss = self.gL2_loss + self.dfL2_loss

        self.g_loss = self.g_loss + self.l2 * self.gL2_loss
        self.d_loss = self.d_loss + self.l2 * self.dfL2_loss

        with tf.variable_scope("g_optimizer", reuse=self.reuse):
            self.g_optimizer = tf.train.AdamOptimizer(self.learning_rate * 5, beta1=0.5)
            self.g_train_step = self.g_optimizer.minimize(self.g_loss, global_step=self.global_step_tensor)

        with tf.variable_scope("d_optimizer", reuse=self.reuse):
            self.d_optimizer = tf.train.AdamOptimizer(self.learning_rate, beta1=0.5)
            self.d_train_step = self.d_optimizer.minimize(self.d_loss, global_step=self.global_step_tensor)

        with tf.variable_scope("ae_loss", reuse=self.reuse):
            self.ae_loss = tf.reduce_mean(self.reconstruction) + self.l2 * self.reconL2_loss  # shape = [None,]

        with tf.variable_scope("ae_optimizer", reuse=self.reuse):
            self.ae_optimizer = tf.train.AdamOptimizer(1.0)
            self.ae_train_step = self.ae_optimizer.minimize(self.ae_loss, global_step=self.global_step_tensor)

    ## ------------------- LOSS: EXPECTED LOWER BOUND ----------------------
    def get_ell(self, x, x_recons):
        """
        Returns the expected log-likelihood of the lower bound.
        For this we use a bernouilli LL.
        """
        # p(x|w)
        return - tf.reduce_sum((x) * tf.log(x_recons + 1e-10) +
                               (1 - x) * tf.log(1 - x_recons + 1e-10), 1)

    '''
    ------------------------------------------------------------------------------
                                     FIT & EVALUATE TENSORS
    ------------------------------------------------------------------------------
    '''

    def partial_fit(self, session, x, z):
        tensors = [self.g_train_step, self.d_train_step, self.g_loss, self.gL2_loss, self.d_loss, self.d_loss_real,
                   self.d_loss_fake, self.dfL2_loss]
        feed_dict = {self.x_batch: x, self.z_batch: z}
        _, _, g_loss, gL2_loss, d_loss, d_loss_real, d_loss_fake, dL2_loss = session.run(tensors, feed_dict=feed_dict)
        return g_loss, gL2_loss, d_loss, d_loss_real, d_loss_fake, dL2_loss

    def evaluate(self, session, x, z):
        tensors = [self.g_loss, self.gL2_loss, self.d_loss, self.d_loss_real, self.d_loss_fake, self.dfL2_loss]
        feed_dict = {self.x_batch: x, self.z_batch: z}
        g_loss, gL2_loss, d_loss, d_loss_real, d_loss_fake, dL2_loss = session.run(tensors, feed_dict=feed_dict)
        return g_loss, gL2_loss, d_loss, d_loss_real, d_loss_fake, dL2_loss

    '''
    ------------------------------------------------------------------------------
                                         GRAPH OPERATIONS
    ------------------------------------------------------------------------------
    '''

    def encode(self, session, inputs):
        tensors = [self.ae_train_step, self.ae_loss]
        feed_dict = {self.x_batch: inputs, self.z_batch: self.z_var}
        _, ae_loss = session.run(tensors, feed_dict=feed_dict)
        ae_loss = np.repeat(np.array([ae_loss]), inputs.shape[0]).reshape((-1, 1))

        return ae_loss

    def decode(self, session, z):
        tensors = [self.x_recons]
        feed_dict = {self.z_batch: z}
        return session.run(tensors, feed_dict=feed_dict)


text:
from math import sqrt, pi
#from mpi_tools import *
import random
import csv
import numpy as np

def random_params(beta, U_range, eps_range, D_range, V_range, N_range, filename="params.csv"):
    U = random.uniform(*U_range)
    eps = random.uniform(*eps_range)
    D = random.uniform(*D_range)
    N = random.randint(*N_range)
    V_list = (V_range[1] - V_range[0]) * np.random.random(N) + V_range[0]
    V_list *= sqrt(2 * D / pi / sum(v**2 for v in V_list))
    #print(2 * D / pi, sum(v**2 for v in V_list))
    e_list = np.sort(np.random.random(N))
    e_list -= sum(v**2 * e for v, e in zip(V_list, e_list)) / (2 * D / pi)
    e_list *= 2 * D / (e_list[-1] - e_list[0])
    #print(0, sum(v**2 * e for v, e in zip(V_list, e_list)))
    #print(2 * D, e_list[-1] - e_list[0])
    params = [U, eps, D, N] + e_list.tolist() + V_list.tolist()
    with open(filename, "ab") as f:
        np.savetxt(f, np.asarray(params)[:,np.newaxis].T, delimiter=",", fmt="%1.4f")

def random_params_cont_bath(beta, U_range, eps_range, D_range, filename="params_cont.csv"):
    U = random.uniform(*U_range)
    eps = random.uniform(*eps_range)
    D = random.uniform(*D_range)
    params = [U, eps, D]
    with open(filename, "ab") as f:
        np.savetxt(f, np.asarray(params)[:,np.newaxis].T, delimiter=",", fmt="%1.4f")

def read_params(filename):
    params_list = []
    with open(filename) as f:
        reader = csv.reader(f)
        for line in reader:
            params = {}
            params["U"] = float(line[0])
            params["eps"] = float(line[1])
            N = int(float(line[3]))
            params["e_list"] = [float(e) for e in line[4:4+N]]
            params["V_list"] = [float(v) for v in line[4+N:4+2*N]]
            params_list.append(params)
    return params_list

def read_params_cont_bath(filename):
    params_list = []
    with open(filename) as f:
        reader = csv.reader(f)
        for line in reader:
            params = {}
            params["U"] = float(line[0])
            params["eps"] = float(line[1])
            params["D"] = float(line[2])
            params_list.append(params)
    return params_list

def name(prefix, beta, n, parent="data/"):
    return parent + prefix + "_beta_" + str(int(beta)) \
            + "_" + str(0 * 10 + n) + ".csv"

beta = 1.

samples_per_file = 1000
files_per_core = 10

if __name__ == "__main__":

    random.seed(get_mpi_rank())
    np.random.seed(get_mpi_rank())
    for n in range(files_per_core):
        for _ in range(samples_per_file):
            random_params(beta=beta,
                          U_range=[1., 8.],
                          eps_range=[-1., 1.],
                          D_range=[2., 8.],
                          V_range=[0.75, 1.25],
                          N_range=[3, 5],
                          filename=name("params", beta, n))

max_stars_repo_name: patryk-kubiczek/learning-GF

max_issues_repo_name: patryk-kubiczek/learning-GF

max_forks_repo_name: patryk-kubiczek/learning-GF


text:
import geopandas as gpd
import pandas as pd
import numpy as np
import movingpandas as mpd
import matplotlib
import matplotlib.pyplot as plt
import seaborn as sns
from tqdm import tqdm
import time
from osm2geojson.helpers import overpass_call
from osm2geojson.main import json2geojson
import warnings
warnings.filterwarnings('ignore')

# based on https://github.com/remisalmon/Strava-local-heatmap-browser
# dark (dark maps, e.g. CartoDB dark_matter), light (light maps, e.g. Stamen Terrain)
# original ("classic" heatmap)

HEATMAP_GRAD = {'dark':{0.0: '#000004',
                        0.1: '#160b39',
                        0.2: '#420a68',
                        0.3: '#6a176e',
                        0.4: '#932667',
                        0.5: '#bc3754',
                        0.6: '#dd513a',
                        0.7: '#f37819',
                        0.8: '#fca50a',
                        0.9: '#f6d746',
                        1.0: '#fcffa4'},
                 'light':{0.0: '#3b4cc0',
                          0.1: '#5977e3',
                          0.2: '#7b9ff9',
                          0.3: '#9ebeff',
                          0.4: '#c0d4f5',
                          0.5: '#dddcdc',
                          0.6: '#f2cbb7',
                          0.7: '#f7ac8e',
                          0.8: '#ee8468',
                          0.9: '#d65244',
                          1.0: '#b40426'},
                 'original':{0: 'black',
                             0.6:'blue',
                             0.7:'lime',
                             0.8:'yellow',
                             0.9:'orange',
                             1:'red'}
                }

def ms_to_km(ms):

    v = (ms) * (60*60)/1000
    return v

def getTopLongestTravelTime(data, top):

    activities = dict()
    for i in range(len(data)):
        activities[i] = [data[i].get_duration(), data[i].tracks[0].type]

    v = sorted(activities.items(), key=lambda item: item[1][0], reverse=True)[:top]
    for j in range(top):
        print("{}° longest activity duration: {}, activity n° {}, type: {}".format(j+1, time.strftime('%H:%M:%S', time.gmtime(v[j][1][0])), v[j][0], v[j][1][1]))

def getTopLongestTravel(data, top):

    activities = dict()
    for i in range(len(data)):
        activities[i] = [data[i].length_3d(), data[i].tracks[0].type]

    v = sorted(activities.items(), key=lambda item: item[1][0], reverse=True)[:top]
    for j in range(top):
        print("{}° longest activity length: {}, activity n° {}, type: {}".format(j+1, round((v[j][1][0]), 2), v[j][0], v[j][1][1]))

def getTopElevationDifference(data, top):

    activities = dict()
    for i in range(len(data)):
        activities[i] = [data[i].get_elevation_extremes().maximum - data[i].get_elevation_extremes().minimum, data[i].tracks[0].type]

    v = sorted(activities.items(), key=lambda item: item[1][0], reverse=True)[:top]
    for j in range(top):
        print("{}° highest elevation difference: {}, activity n° {}, type: {}".format(j+1, round((v[j][1][0]), 2), v[j][0], v[j][1][1]))

def toList(activities) -> list:

    target = []
    for activity in tqdm(activities):
        #print(activity)
        data = []
        for point_idx, point in enumerate(activity.tracks[0].segments[0].points):
            data.append([point.longitude, point.latitude, point.elevation, point.time])

        cols = ['longitude', 'latitude', 'elevation', 'time']
        gpx_track = pd.DataFrame(data, columns=cols)
        gpx_track.time = gpx_track.time.apply(lambda x: x.replace(tzinfo=None))
        target.append(gpx_track)

    return target

def toGdfList(activities_dfList) -> list:

    target_gdfList = []
    for i in tqdm(range(len(activities_dfList))):
        geo_df = gpd.GeoDataFrame(activities_dfList[i], crs=4326, geometry=gpd.points_from_xy(activities_dfList[i].longitude, activities_dfList[i].latitude, activities_dfList[i].elevation))
        # geo_df.set_index('time', drop=True, inplace=True)
        target_gdfList.append(geo_df)

    return target_gdfList

def getTrajList(geo_dfList) -> list:

    trajectories = []
    for i in range(len(geo_dfList)):
        trajectory = mpd.Trajectory(geo_dfList[i], i)
        trajectories.append(trajectory)

    return trajectories

def get_boundary(city: str):

    boundary = gpd.GeoDataFrame()
    print("Getting boundaries for {}".format(city))

    query = '''
        [out:json];
        area[name="{}"][admin_level=8][boundary=administrative]->.target;
        area["name"="Italia"][boundary=administrative]->.wrap;
        rel(pivot.target)(area.wrap);

        out geom;
    '''.format(city)

    ti = time.time()
    result = overpass_call(query)
    res = json2geojson(result)
    boundary = boundary.append(gpd.GeoDataFrame.from_features(res), ignore_index=True)
    boundary.set_crs(epsg=4326, inplace=True)
    #coords = [i for i in boundary.to_crs(epsg=3587).geometry[0].envelope.exterior.coords]

    if len(boundary) <= 0:
        raise ConnectionError("No boundaries were found for {}. Try with another city or check your Overpass query limit.".format(city))
    else:
        print("Extracted boundaries for {}. Time elapsed: {} s".format(city, round(time.time()-ti, 2)))

    #return min(coords[0][0], coords[2][0]), max(coords[0][0], coords[2][0]), min(coords[0][1], coords[2][1]), max(coords[0][1], coords[2][1])
    return boundary

def getStopElevationDiff(stop_points, start) -> gpd.GeoDataFrame:

    stop_points['elevation_diff'] = 0.0
    stop_points['time_diff'] = stop_points.end_time[0] - stop_points.start_time[0]

    for i in range(len(stop_points)):
        if i==0:
            stop_points.elevation_diff[i] = np.round(stop_points.geometry[i].z - start.df.elevation[0], 2)
            stop_points.time_diff[i] = stop_points.start_time[i] - start.df.index[0]
        else:
            stop_points.elevation_diff[i] = np.round(stop_points.geometry[i].z - stop_points.geometry[i-1].z, 2)
            stop_points.time_diff[i] = stop_points.start_time[i] - stop_points.end_time[i-1]

    return stop_points

def plotHex(data, data_extent, basemap, basemap_extent, hex):

    f, ax = plt.subplots(1, figsize=(12, 12))
    ax.imshow(basemap, extent=basemap_extent, interpolation='bilinear')
    # Generate and add hexbin with hex hexagons in each
    # dimension, no borderlines, 60% transparency (alpha),
    # and the reverse Oranges colormap
    hb = ax.hexbin(x=data['x'],
                y=data['y'],
                gridsize=hex, linewidths=0,
                alpha=0.6, cmap='Reds')
    ax.axis(data_extent)
    plt.colorbar(hb)
    ax.set_axis_off()
    plt.show()

def plotKDE(data, data_extent, basemap, basemap_extent, grad):

    f, ax = plt.subplots(1, figsize=(15, 15))
    ax.imshow(basemap, extent=basemap_extent, interpolation='bilinear')
    # Generate and add KDE with a shading of grad gradients
    # coloured contours, 30% of transparency (alpha),
    # and the reverse Reds colormap
    sns.kdeplot(x=data['x'], y=data['y'],
                    n_levels=grad, shade=True,
                    alpha=0.3, cmap='Reds')
    ax.axis(data_extent)
    ax.set_axis_off()
    plt.show()

def plotClusters(data, labels, basemap, basemap_extent):

    colors = ['grey','yellow','blue','purple', 'red', 'orange', 'green']
    f, ax = plt.subplots(1, figsize=(12, 12))
    ax.imshow(basemap, extent=basemap_extent, interpolation='bilinear')
    ax.scatter(data['x'], data['y'], \
    c=labels, cmap=matplotlib.colors.ListedColormap(colors), linewidth=0)

    ax.set_axis_off()
    plt.show()

def plotRunComparison(data):

    fig, axs = plt.subplots(3, sharex=True, figsize=(15,10))
    fig.suptitle('Running performances comparison')
    x = data.date
    axs[0].plot(x, data.length, color='orange')
    axs[0].set_ylabel('km', color='black')
    axs[0].grid(b=True, which='both', axis='both')
    axs[1].plot(x, data.avgPace, color='blue')
    axs[1].set_ylabel('km/h', color='black')
    axs[1].grid(b=True, which='both', axis='both')
    axs[2].plot(x, data.vo2MaxValue, color='red')
    axs[2].set_ylabel('ml/kg/min', color='black')
    axs[2].grid(b=True, which='both', axis='both')
    axs[2].set_xlabel('date')

    for ax in axs:
        ax.xaxis.set_major_locator(plt.MaxNLocator(25))
        ax.label_outer()

    plt.xticks(rotation='vertical')

    plt.show()