text:
import cv2
import numpy as np

# TODO: detect windows
# TODO: detect doors
# Calculate (actual) size of appartment
# TODO: text detection

"""
Detect
This file contains functions used when detecting and calculating shapes in images.

FloorplanToBlender3d
Copyright (C) 2019 Daniel Westberg
"""

def wall_filter(gray):
    """
    Filter walls
    Filter out walls from a grayscale image
    @Param image
    @Return image of walls
    """
    ret, thresh = cv2.threshold(gray,0,255,cv2.THRESH_BINARY_INV+cv2.THRESH_OTSU)

    # noise removal
    kernel = np.ones((3,3),np.uint8)
    opening = cv2.morphologyEx(thresh,cv2.MORPH_OPEN,kernel, iterations = 2)

    sure_bg = cv2.dilate(opening,kernel,iterations=3)

    dist_transform = cv2.distanceTransform(opening,cv2.DIST_L2,5)
    ret, sure_fg = cv2.threshold(0.5*dist_transform,0.2*dist_transform.max(),255,0)

    sure_fg = np.uint8(sure_fg)
    unknown = cv2.subtract(sure_bg,sure_fg)

    return unknown

def detectPreciseBoxes(detect_img, output_img = None, color = [100,100,0]):
    """
    Detect corners with boxes in image with high precision
    @Param detect_img image to detect from @mandatory
    @Param output_img image for output
    @Param color to set on output
    @Return corners(list of boxes), output image
    @source https://stackoverflow.com/questions/50930033/drawing-lines-and-distance-to-them-on-image-opencv-python
    """
    res = []

    contours, hierarchy = cv2.findContours(detect_img,cv2.RETR_EXTERNAL,cv2.CHAIN_APPROX_SIMPLE)
    #area = sorted(contours, key=cv2.contourArea, reverse=True)

    largest_contour_area = 0
    for cnt in contours:
        largest_contour_area = cv2.contourArea(cnt)
        largest_contour = cnt

        epsilon = 0.001*cv2.arcLength(largest_contour,True)
        approx = cv2.approxPolyDP(largest_contour,epsilon,True)
        if output_img is not None:
            final = cv2.drawContours(output_img, [approx], 0, color)
        res.append(approx)

    return res, output_img

def remove_noise(img, noise_removal_threshold):
    """
    Remove noise from image and return mask
    Help function for finding room
    @Param img @mandatory image to remove noise from
    @Param noise_removal_threshold @mandatory threshold for noise
    @Return return new mask of image
    """
    img[img < 128] = 0
    img[img > 128] = 255
    contours, _ = cv2.findContours(~img, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
    mask = np.zeros_like(img)
    for contour in contours:
        area = cv2.contourArea(contour)
        if area > noise_removal_threshold:
            cv2.fillPoly(mask, [contour], 255)
    return mask

def find_corners_and_draw_lines(img, corners_threshold, room_closing_max_length):
    """
    Finds corners and draw lines from them
    Help function for finding room
    @Param image input image
    @Param corners_threshold threshold for corner distance
    @Param room_closing_max_length threshold for room max size
    @Return output image
    """
    # Detect corners (you can play with the parameters here)
    kernel = np.ones((1,1),np.uint8)
    dst = cv2.cornerHarris(img ,2,3,0.04)
    dst = cv2.erode(dst,kernel, iterations = 10)
    corners = dst > corners_threshold * dst.max()

    # Draw lines to close the rooms off by adding a line between corners on the same x or y coordinate
    # This gets some false positives.
    # You could try to disallow drawing through other existing lines for example.
    for y,row in enumerate(corners):
        x_same_y = np.argwhere(row)
        for x1, x2 in zip(x_same_y[:-1], x_same_y[1:]):

            if x2[0] - x1[0] < room_closing_max_length:
                color = 0
                cv2.line(img, (x1[0], y), (x2[0], y), color, 1)

    for x,col in enumerate(corners.T):
        y_same_x = np.argwhere(col)
        for y1, y2 in zip(y_same_x[:-1], y_same_x[1:]):
            if y2[0] - y1[0] < room_closing_max_length:
                color = 0
                cv2.line(img, (x, y1[0]), (x, y2[0]), color, 1)
    return img

def mark_outside_black(img, mask):
    """
    Mark white background as black
    @Param @mandatory img image input
    @Param @mandatory mask mask to use
    @Return image, mask
    """
    # Mark the outside of the house as black
    contours, _ = cv2.findContours(~img, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
    contour_sizes = [(cv2.contourArea(contour), contour) for contour in contours]
    biggest_contour = max(contour_sizes, key=lambda x: x[0])[1]
    mask = np.zeros_like(mask)
    cv2.fillPoly(mask, [biggest_contour], 255)
    img[mask == 0] = 0
    return img, mask

def find_rooms(img, noise_removal_threshold=50, corners_threshold=0.01,
               room_closing_max_length=130,
               gap_in_wall_min_threshold=5000):

    """
    I have copied and changed this function some...

    origin from
    https://stackoverflow.com/questions/54274610/crop-each-of-them-using-opencv-python

    @param img: grey scale image of rooms, already eroded and doors removed etc.
    @param noise_removal_threshold: Minimal area of blobs to be kept.
    @param corners_threshold: Threshold to allow corners. Higher removes more of the house.
    @param room_closing_max_length: Maximum line length to add to close off open doors.
    @param gap_in_wall_threshold: Minimum number of pixels to identify component as room instead of hole in the wall.
    @return: rooms: list of numpy arrays containing boolean masks for each detected room
             colored_house: A colored version of the input image, where each room has a random color.
    """
    assert 0 <= corners_threshold <= 1
    # Remove noise left from door removal

    mask = remove_noise(img, noise_removal_threshold)
    img = ~mask

    find_corners_and_draw_lines(img,corners_threshold,room_closing_max_length)

    img, mask = mark_outside_black(img, mask)

    # Find the connected components in the house
    ret, labels = cv2.connectedComponents(img)
    img = cv2.cvtColor(img,cv2.COLOR_GRAY2RGB)
    unique = np.unique(labels)
    rooms = []
    for label in unique:
        component = labels == label
        if img[component].sum() == 0 or np.count_nonzero(component) < gap_in_wall_min_threshold:
            color = 0
        else:
            rooms.append(component)
            color = np.random.randint(0, 255, size=3)
        img[component] = color
    return rooms, img

def detectAndRemovePreciseBoxes(detect_img, output_img = None, color = [255, 255, 255]):
    """
    Remove contours of detected walls from image
    @Param detect_img image to detect from @mandatory
    @Param output_img image for output
    @Param color to set on output
    @Return list of boxes, actual image
    @source https://stackoverflow.com/questions/50930033/drawing-lines-and-distance-to-them-on-image-opencv-python
    """

    res = []

    contours, hierarchy = cv2.findContours(detect_img,cv2.RETR_EXTERNAL,cv2.CHAIN_APPROX_SIMPLE)
    #area = sorted(contours, key=cv2.contourArea, reverse=True)

    largest_contour_area = 0
    for cnt in contours:
        largest_contour_area = cv2.contourArea(cnt)
        largest_contour = cnt

        epsilon = 0.001*cv2.arcLength(largest_contour,True)
        approx = cv2.approxPolyDP(largest_contour,epsilon,True)
        if output_img is not None:
            cv2.drawContours( output_img,  [approx], -1, color, -1);
        res.append(approx)

    return res, output_img

def detectOuterContours(detect_img, output_img = None, color = [255, 255, 255]):
    """
    Get the outer side of floorplan, used to get ground
    @Param detect_img image to detect from @mandatory
    @Param output_img image for output
    @Param color to set on output
    @Return approx, box
    @Source https://stackoverflow.com/questions/50930033/drawing-lines-and-distance-to-them-on-image-opencv-python
    """
    ret, thresh = cv2.threshold(detect_img, 230, 255, cv2.THRESH_BINARY_INV)

    contours, hierarchy = cv2.findContours(thresh.copy(), cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)

    largest_contour_area = 0
    for cnt in contours:
        if (cv2.contourArea(cnt) > largest_contour_area):
            largest_contour_area = cv2.contourArea(cnt)
            largest_contour = cnt

    epsilon = 0.001*cv2.arcLength(largest_contour,True)
    approx = cv2.approxPolyDP(largest_contour,epsilon,True)
    if output_img is not None:
        final = cv2.drawContours(output_img, [approx], 0, color)
    return approx, output_img

def rectContains(rect,pt):
    """
    Count if Rect contains point
    @Param rect rectangle
    @Param pt point
    @Return boolean
    @source: https://stackoverflow.com/questions/33065834/how-to-detect-if-a-point-is-contained-within-a-bounding-rect-opecv-python
    """
    return rect[0] < pt[0] < rect[0]+rect[2] and rect[1] < pt[1] < rect[1]+rect[3]

'''
Currently none used code below here!, outcommented to avoid confusion with contributors.
'''

def find_details(img, noise_removal_threshold=50, corners_threshold=0.01,
               room_closing_max_length=130, gap_in_wall_max_threshold=5000,
               gap_in_wall_min_threshold=10):

    """
    !!! Currently not used in IMPLEMENTATION !!!
    I have copied and changed this function some...

    origin from
    https://stackoverflow.com/questions/54274610/crop-each-of-them-using-opencv-python

    @Param img: grey scale image of rooms, already eroded and doors removed etc.
    @Param noise_removal_threshold: Minimal area of blobs to be kept.
    @Param corners_threshold: Threshold to allow corners. Higher removes more of the house.
    @Param room_closing_max_length: Maximum line length to add to close off open doors.
    @Param gap_in_wall_threshold: Minimum number of pixels to identify component as room instead of hole in the wall.
    @Return: rooms: list of numpy arrays containing boolean masks for each detected room
             colored_house: A colored version of the input image, where each room has a random color.
    """
    assert 0 <= corners_threshold <= 1
    # Remove noise left from door removal

    mask = remove_noise(img, noise_removal_threshold)
    img = ~mask

    find_corners_and_draw_lines(img,corners_threshold,room_closing_max_length)

    img, mask = mark_outside_black(img, mask)

    # Find the connected components in the house
    ret, labels = cv2.connectedComponents(img)
    img = cv2.cvtColor(img,cv2.COLOR_GRAY2RGB)
    unique = np.unique(labels)
    details = []
    for label in unique:
        component = labels == label
        if img[component].sum() == 0 or np.count_nonzero(component) < gap_in_wall_min_threshold or np.count_nonzero(component) > gap_in_wall_max_threshold:
            color = 0
        else:
            details.append(component)
            color = np.random.randint(0, 255, size=3)

        img[component] = color

    return details, img

max_stars_repo_name: rbg-research/Floor-Plan-Detection

max_issues_repo_name: rbg-research/Floor-Plan-Detection

max_forks_repo_name: rbg-research/Floor-Plan-Detection


text:
import os, stat, time
from pathlib import Path

import numpy as np

from ..system import System

def read_fifo(fin, dtype, count):
    buffer = fin.read(int(dtype[-1]) * count)
    return np.frombuffer(buffer, dtype=dtype, count=count)

class IOFifo(object):
    def __init__(self, cwd=None):
        self.mode = "fifo"
        self.cwd = cwd

    def load_system(self, input, system=None, step=None):

        self.cwd = self.cwd or Path(input).parent

        if step is None:
            step = 0

        self._step = np.asarray(step)

        assert stat.S_ISFIFO(os.stat(input).st_mode)

        self._fin = open(input, "rb")

        self._n_atoms, self._n_qm_atoms, qm_charge, qm_mult, self._pbc = read_fifo(self._fin, dtype="i4", count=5)
        self._system = system or System(self._n_atoms, self._n_qm_atoms, qm_charge=qm_charge, qm_mult=qm_mult)

        return self._system

    def return_results(self, energy, forces, output=None):
        assert self._system is not None

        output = output or Path(self._fin.name).with_suffix('.out')

        try:
            os.mkfifo(output)
        except:
            pass

        self._system.atoms.charges[:] = read_fifo(self._fin, dtype="f8", count=self._n_atoms)
        self._system.qm.atoms.elements[:] = read_fifo(self._fin, dtype="i4", count=self._n_qm_atoms)

        if self._pbc > 0:
            cell_basis = read_fifo(self._fin, dtype="f8", count=9).reshape(3, 3).copy()
            cell_basis[np.isclose(cell_basis, 0.0)] = 0.0
            self._system.cell_basis[:] = cell_basis

        # First cycle
        self._step[()] = read_fifo(self._fin, dtype="i4", count=1)

        if self._pbc == 2 :
            cell_basis = read_fifo(self._fin, dtype="f8", count=9).reshape(3, 3).copy()
            cell_basis[np.isclose(cell_basis, 0.0)] = 0.0
            self._system.cell_basis[:] = cell_basis

        self._system.atoms.positions[:] = read_fifo(self._fin, dtype="f8", count=self._n_atoms * 3).reshape(3, self._n_atoms)

        self._fout = open(output, "wb")
        self._fout.write(energy.tobytes())
        self._fout.write(forces.tobytes(order="F"))

        while True:

            try:
                _step = read_fifo(self._fin, dtype="i4", count=1)
                self._step[()] = _step
            except:
                break

            if self._pbc == 2 :
                cell_basis = read_fifo(self._fin, dtype="f8", count=9).reshape(3, 3).copy()
                cell_basis[np.isclose(cell_basis, 0.0)] = 0.0
                self._system.cell_basis[:] = cell_basis

            self._system.atoms.positions[:] = read_fifo(self._fin, dtype="f8", count=self._n_atoms * 3).reshape(3, self._n_atoms)

            self._fout.write(energy.tobytes())
            self._fout.write(forces.tobytes(order="F"))

    @staticmethod
    def save_input(input):
        """Preserve the input file passed from the driver."""
        pass

ext: py

max_stars_repo_name: panxl/qmhub

max_issues_repo_name: panxl/qmhub

max_forks_repo_name: panxl/qmhub

reason: import numpy


text:
#!/usr/bin/env python
# -*- coding: utf8 -*-
# *****************************************************************
# **       PTS -- Python Toolkit for working with SKIRT          **
# **       © Astronomical Observatory, Ghent University          **
# *****************************************************************

## \package pts.modeling.core.sed Contains the SED class.

# -----------------------------------------------------------------

# Ensure Python 3 compatibility
from __future__ import absolute_import, division, print_function

# Import standard modules
import numpy as np

# Import astronomical modules
from astropy.units import Unit, spectral
from astropy.table import Table

# Import the relevant PTS classes and modules
from ...core.tools import tables
from ...core.tools import introspection
from ...core.tools import filesystem as fs
from ...core.basics.errorbar import ErrorBar
from ...core.basics.filter import Filter

# -----------------------------------------------------------------

class IntrinsicSED(object):

    """
    This class ...
    """

    def __init__(self):

        """
        This function ...
        """

        # Attributes
        self.table = None

    # -----------------------------------------------------------------

    def wavelengths(self, unit=None, asarray=False, add_unit=True):

        """
        This function ...
        :param unit:
        :param asarray:
        :param add_unit:
        :return:
        """

        if asarray: return tables.column_as_array(self.table["Wavelength"], unit=unit)
        else: return tables.column_as_list(self.table["Wavelength"], unit=unit, add_unit=add_unit)

    # -----------------------------------------------------------------

    def luminosities(self, unit=None, asarray=False, add_unit=True):

        """
        This function ...
        :param unit:
        :param asarray:
        :param add_unit:
        :return:
        """

        if asarray: return tables.column_as_array(self.table["Luminosity"], unit=unit)
        else: return tables.column_as_list(self.table["Luminosity"], unit=unit, add_unit=add_unit)

    # -----------------------------------------------------------------

    @classmethod
    def from_file(cls, path, skiprows=0):

        """
        This function ...
        :param path:
        :param skiprows:
        :return:
        """

        # Create a new SED
        sed = cls()

        wavelength_column, luminosity_column = np.loadtxt(path, dtype=float, unpack=True, skiprows=skiprows)
        sed.table = tables.new([wavelength_column, luminosity_column], ["Wavelength", "Luminosity"])
        sed.table["Wavelength"].unit = Unit("micron")
        sed.table["Luminosity"].unit = Unit("W/micron")

        # Return the SED
        return sed

    # -----------------------------------------------------------------

    @classmethod
    def from_luminosities(cls, wavelengths, luminosities, wavelength_unit="micron", luminosity_unit="W/micron"):

        """
        This function ...
        :return:
        """

        # Create a new SED
        sed = cls()

        sed.table = tables.new([wavelengths, luminosities], ["Wavelength", "Luminosity"])
        sed.table["Wavelength"].unit = wavelength_unit
        sed.table["Luminosity"].unit = luminosity_unit

        # Return the SED
        return sed

# -----------------------------------------------------------------

class ObservedSED(object):

    """
    This class ...
    """

    def __init__(self):

        """
        This function ...
        :return:
        """

        # Attributes
        self.table = Table(names=["Observatory", "Instrument", "Band", "Wavelength", "Flux", "Error-", "Error+"],
                           dtype=('S10', 'S10', 'S10', 'f8', 'f8', 'f8', 'f8'))
        self.table["Wavelength"].unit = Unit("micron")
        self.table["Flux"].unit = Unit("Jy")
        self.table["Error-"].unit = Unit("Jy")
        self.table["Error+"].unit = Unit("Jy")

    # -----------------------------------------------------------------

    @classmethod
    def from_file(cls, path):

        """
        This function ...
        :param path:
        :return:
        """

        # Create a new observed SED
        sed = cls()

        #names = ["Observatory", "Instrument", "Band", "Wavelength", "Flux", "Error-", "Error+"]
        #observatory_column, instrument_column, band_column, wavelength_column, flux_column, error_min_column, error_plus_column = np.loadtxt(path, unpack=True, dtype=str)
        #wavelength_column = wavelength_column.astype(float)
        #flux_column = flux_column.astype(float)
        #error_min_column = error_min_column.astype(float)
        #error_plus_column = error_plus_column.astype(float)
        #sed.table = tables.new([observatory_column, instrument_column, band_column, wavelength_column, flux_column, error_min_column, error_plus_column], names)
        #sed.table["Wavelength"].unit = "micron"
        #sed.table["Flux"].unit = "Jy"
        #sed.table["Error-"].unit = "Jy"
        #sed.table["Error+"].unit = "Jy"

        # New
        sed.table = tables.from_file(path, format="ascii.ecsv")

        # Return the observed SED
        return sed

    # -----------------------------------------------------------------

    @classmethod
    def from_caapr(cls, path):

        """
        This function ...
        :return:
        """

        # Create a new observed SED
        sed = cls()

        # Load the table
        caapr_table = tables.from_file(path, format="csv")

        fluxes = dict()
        errors = dict()

        # Loop over the columns of the table
        for colname in caapr_table.colnames:

            if colname == "name": continue
            if "ERR" in colname:
                instrument_band = colname.split("_ERR")[0]
                error = abs(caapr_table[colname][0])
                if not np.isnan(error): errors[instrument_band] = error
            else:
                instrument_band = colname
                flux = caapr_table[colname][0]
                if not np.isnan(flux): fluxes[instrument_band] = flux

        observatory_column = []
        instrument_column = []
        band_column = []
        wavelength_column = []
        flux_column = []
        fluxerrmin_column = []
        fluxerrmax_column = []

        for instrument_band in fluxes:

            if not instrument_band in errors: raise ValueError("No error for " + instrument_band)

            flux = fluxes[instrument_band]
            error = errors[instrument_band]

            instrument = instrument_band.split("_")[0]
            band = instrument_band.split("_")[1]

            # Create filter
            fltr = Filter.from_string(instrument + " " + band)

            # Get filter properties
            observatory = fltr.observatory
            instrument = fltr.instrument
            band = fltr.band
            wavelength = fltr.pivotwavelength()

            # Add entry to the columns
            observatory_column.append(observatory)
            instrument_column.append(instrument)
            band_column.append(band)
            wavelength_column.append(wavelength)
            flux_column.append(flux)
            fluxerrmin_column.append(-error)
            fluxerrmax_column.append(error)

        # Create the SED table
        data = [observatory_column, instrument_column, band_column, wavelength_column, flux_column, fluxerrmin_column, fluxerrmax_column]
        names = ["Observatory", "Instrument", "Band", "Wavelength", "Flux", "Error-", "Error+"]
        sed.table = tables.new(data, names)
        sed.table["Wavelength"].unit = "micron"
        sed.table["Flux"].unit = "Jy"
        sed.table["Error-"].unit = "Jy"
        sed.table["Error+"].unit = "Jy"

        # Return the SED
        return sed

    # -----------------------------------------------------------------

    def add_entry(self, fltr, flux, error):

        """
        This function ...
        :param fltr:
        :param flux:
        :param error:
        :return:
        """

        self.table.add_row([fltr.observatory, fltr.instrument, fltr.band, fltr.pivotwavelength(), flux, error.lower, error.upper])

    # -----------------------------------------------------------------

    def instruments(self):

        """
        This function ...
        :return:
        """

        return tables.column_as_list(self.table["Instrument"])

    # -----------------------------------------------------------------

    def bands(self):

        """
        This function ...
        :return:
        """

        return tables.column_as_list(self.table["Band"])

    # -----------------------------------------------------------------

    def filters(self):

        """
        This function ...
        :return:
        """

        # Initialize
        filters = []

        # Loop over all entries
        for i in range(len(self.table)):

            # Get the instrument and band
            instrument = self.table["Instrument"][i]
            band = self.table["Band"][i]

            # Create the filter
            fltr = Filter.from_instrument_and_band(instrument, band)

            # Add the filter to the list
            filters.append(fltr)

        # Return the list of filters
        return filters

    # -----------------------------------------------------------------

    def wavelengths(self, unit=None, asarray=False, add_unit=True):

        """
        This function ...
        :param unit:
        :param asarray:
        :param add_unit:
        :return:
        """

        if asarray: return tables.column_as_array(self.table["Wavelength"], unit=unit)
        else: return tables.column_as_list(self.table["Wavelength"], unit=unit, add_unit=add_unit)

    # -----------------------------------------------------------------

    def fluxes(self, unit=None, asarray=False, add_unit=True):

        """
        This function ...
        :param unit:
        :param asarray:
        :param add_unit:
        :return:
        """

        if asarray: return tables.column_as_array(self.table["Flux"], unit=unit)
        else: return tables.column_as_list(self.table["Flux"], unit=unit, add_unit=add_unit)

    # -----------------------------------------------------------------

    def errors(self, unit=None, add_unit=True):

        """
        This function ...
        :param unit:
        :param add_unit:
        :return:
        """

        return tables.columns_as_objects([self.table["Error-"], self.table["Error+"]], ErrorBar, unit=unit, add_unit=add_unit)

    # -----------------------------------------------------------------

    def flux_for_filter(self, fltr, unit=None, add_unit=True):

        """
        This function ...
        :param fltr:
        :param unit:
        :param add_unit:
        :return:
        """

        return self.flux_for_band(fltr.instrument, fltr.band, unit, add_unit)

    # -----------------------------------------------------------------

    def flux_for_band(self, instrument, band, unit=None, add_unit=True):

        """
        This function ...
        :param instrument:
        :param band:
        :param unit:
        :param add_unit:
        :return:
        """

        has_unit = self.table["Flux"].unit is not None
        has_mask = hasattr(self.table["Flux"], "mask")

        # If unit has to be converted, check whether the original unit is specified
        if not has_unit and unit is not None: raise ValueError(
            "Cannot determine the unit of the flux column so values cannot be converted to " + str(unit))

        # Loop over all the entries in the table
        for i in range(len(self.table)):

            instrument_entry = self.table["Instrument"][i]
            band_entry = self.table["Band"][i]

            if not (instrument_entry == instrument and band_entry == band): continue

            if has_unit:

                # Add the unit initially to be able to convert
                flux = self.table["Flux"][i] * self.table["Flux"].unit

                # If a target unit is specified, convert
                if unit is not None: flux = flux.to(unit).value * Unit(unit)

                if not add_unit: flux = flux.value

            else: flux = self.table["Flux"][i]

            return flux

        # If no match is found, return None
        return None

    # -----------------------------------------------------------------

    def error_for_filter(self, fltr, unit=None, add_unit=True):

        """
        This function ...
        :param fltr:
        :param unit:
        :param add_unit:
        :return:
        """

        return self.error_for_band(fltr.instrument, fltr.band, unit, add_unit)

    # -----------------------------------------------------------------

    def error_for_band(self, instrument, band, unit=None, add_unit=True):

        """
        This function ...
        :param instrument:
        :param band:
        :param unit:
        :param add_unit:
        :return:
        """

        has_unit = self.table["Error-"].unit is not None and self.table["Error+"].unit is not None
        has_mask = hasattr(self.table["Error-"], "mask")
        assert has_mask == hasattr(self.table["Error+"], "mask")

        # If unit has to be converted, check whether the original unit is specified
        if not has_unit and unit is not None: raise ValueError(
            "Cannot determine the unit of the error columns so values cannot be converted to " + str(unit))

        # Loop over all the entries in the table
        for i in range(len(self.table)):

            instrument_entry = self.table["Instrument"][i]
            band_entry = self.table["Band"][i]

            if not (instrument_entry == instrument and band_entry == band): continue

            if has_unit:

                # Add the unit initially to be able to convert
                error_min = self.table["Error-"][i] * self.table["Error-"].unit
                error_plus = self.table["Error+"][i] * self.table["Error+"].unit

                # If a target unit is specified, convert
                if unit is not None:

                    error_min = error_min.to(unit).value * Unit(unit)
                    error_plus = error_plus.to(unit).value * Unit(unit)

                if not add_unit:

                    error_min = error_min.value
                    error_plus = error_plus.value

                error = ErrorBar(error_min, error_plus)

            else: error = ErrorBar(self.table["Error-"][i], self.table["Error+"][i])

            return error

        # If no match is found, return None
        return None

    # -----------------------------------------------------------------

    def save(self, path):

        """
        This function ...
        :param path:
        :return:
        """

        # Sort the table by wavelength
        self.table.sort("Wavelength")

        # Write the observed SED
        tables.write(self.table, path, format="ascii.ecsv")

# -----------------------------------------------------------------

class SED(object):

    """
    This class...
    """

    def __init__(self, wavelength_unit="micron", flux_unit="Jy"):

        """
        The constructor ...
        :return:
        """

        # Attributes
        self.table = Table(names=["Wavelength", "Flux", "Error-", "Error+"], dtype=('f8', 'f8', 'f8', 'f8'))
        self.table["Wavelength"].unit = Unit(wavelength_unit)
        self.table["Flux"].unit = Unit(flux_unit)
        self.table["Error-"].unit = Unit(flux_unit)
        self.table["Error+"].unit = Unit(flux_unit)

    # -----------------------------------------------------------------

    def add_entry(self, wavelength, flux, error=None):

        """
        This function ...
        :param wavelength:
        :param flux:
        :param error:
        :return:
        """

        wavelength_unit = self.table["Wavelength"].unit
        wavelength = wavelength.to(wavelength_unit).value

        flux_unit = self.table["Flux"].unit
        flux = flux.to(flux_unit).value

        error_lower = error.lower.to(flux_unit).value if error is not None else None
        error_upper = error.upper.to(flux_unit).value if error is not None else None

        self.table.add_row([wavelength, flux, error_lower, error_upper])

    # -----------------------------------------------------------------

    def wavelengths(self, unit=None, asarray=False, add_unit=True):

        """
        This function ...
        :param unit:
        :param asarray:
        :param add_unit:
        :return:
        """

        if asarray: return tables.column_as_array(self.table["Wavelength"], unit=unit)
        else: return tables.column_as_list(self.table["Wavelength"], unit=unit, add_unit=add_unit)

    # -----------------------------------------------------------------

    def fluxes(self, unit=None, asarray=False, add_unit=True):

        """
        This function ...
        :param unit:
        :param asarray:
        :param add_unit:
        :return:
        """

        if asarray: return tables.column_as_array(self.table["Flux"], unit=unit)
        else: return tables.column_as_list(self.table["Flux"], unit=unit, add_unit=add_unit)

    # -----------------------------------------------------------------

    def errors(self, unit=None, add_unit=True):

        """
        This function ...
        :param unit:
        :param add_unit:
        :return:
        """

        return tables.columns_as_objects([self.table["Error-"], self.table["Error+"]], ErrorBar, unit=unit, add_unit=add_unit)

    # -----------------------------------------------------------------

    @property
    def has_errors(self):

        """
        This function ...
        :return:
        """

        return "Error-" in self.table.colnames and "Error+" in self.table.colnames

    # -----------------------------------------------------------------

    @classmethod
    def from_file(cls, path):

        """
        This function ...
        :param path:
        :return:
        """

        # Create a new SED
        sed = cls()

        # New
        sed.table = tables.from_file(path, format="ascii.ecsv")

        # Return the SED
        return sed

    # -----------------------------------------------------------------

    @classmethod
    def from_skirt(cls, path, skiprows=0, contribution="total"):

        """
        This function ...
        :param path:
        :param skiprows:
        :param contribution:
        :return:
        """

        # Create a new SED
        sed = cls()

        # SEDInstrument:
        # column 1: lambda (micron)
        # column 2: total flux; lambda*F_lambda (W/m2)

        # From FullInstrument:
        # column 1: lambda (micron)
        # column 2: total flux; lambda*F_lambda (W/m2)
        # column 3: direct stellar flux; lambda*F_lambda (W/m2)
        # column 4: scattered stellar flux; lambda*F_lambda (W/m2)
        # column 5: total dust emission flux; lambda*F_lambda (W/m2)
        # column 6: dust emission scattered flux; lambda*F_lambda (W/m2)
        # column 7: transparent flux; lambda*F_lambda (W/m2)

        from ..preparation import unitconversion

        # Open the SED table
        #sed.table = tables.from_file(path, format="ascii.no_header") # sometimes doesn't work ?? why ??
        #sed.table.rename_column("col1", "Wavelength")
        #sed.table.rename_column("col2", "Flux")

        if contribution == "total": columns = (0,1)
        elif contribution == "direct": columns = (0,2)
        elif contribution == "scattered": columns = (0,3)
        elif contribution == "dust": columns = (0,4)
        elif contribution == "dustscattered": columns = (0,5)
        elif contribution == "transparent": columns = (0,6)
        else: raise ValueError("Wrong value for 'contribution': should be 'total', 'direct', 'scattered', 'dust', 'dustscattered' or 'transparent'")

        wavelength_column, flux_column = np.loadtxt(path, dtype=float, unpack=True, skiprows=skiprows, usecols=columns)

        sed.table = tables.new([wavelength_column, flux_column], ["Wavelength", "Flux"])
        sed.table["Wavelength"].unit = Unit("micron")

        jansky_column = []

        for i in range(len(sed.table)):

            # Get the flux density in W / m2 and the wavelength in micron
            neutral_fluxdensity = sed.table["Flux"][i] * Unit("W/m2")
            wavelength = sed.table["Wavelength"][i] * Unit("micron")

            # Convert to Jansky (2 methods give same result)
            #jansky_ = unitconversion.neutral_fluxdensity_to_jansky(neutral_fluxdensity, wavelength)
            jansky = (neutral_fluxdensity / wavelength.to("Hz", equivalencies=spectral())).to("Jy").value

            # Add the fluxdensity in Jansky to the new column
            jansky_column.append(jansky)

        # Add the flux column in Jansky
        sed.table.remove_column("Flux")
        sed.table["Flux"] = jansky_column
        sed.table["Flux"].unit = "Jy"

        # Return the SED
        return sed

    # -----------------------------------------------------------------

    def save(self, path):

        """
        This function ...
        :param path
        :return:
        """

        # Write the SED table to file
        tables.write(self.table, path, format="ascii.ecsv")

# -----------------------------------------------------------------

seds_path = fs.join(introspection.pts_dat_dir("modeling"), "seds")

# -----------------------------------------------------------------

def load_example_mappings_sed():

    """
    This function ...
    :return:
    """

    # Determine the path to the SED file
    sed_path = fs.join(seds_path, "mapsed.dat")

    # Get the data
    wavelength_column, flux_column = np.loadtxt(sed_path, usecols=(0, 1), unpack=True)

    # Create an SED instance
    sed = SED()

    # Set the columns
    sed.table = tables.new([wavelength_column, flux_column], ["Wavelength", "Flux"])
    sed.table["Wavelength"].unit = "micron"
    sed.table["Flux"].unit = "W/m2" # = lambda * F_Lambda !

    # Return the SED
    return sed

# -----------------------------------------------------------------

def load_example_bruzualcharlot_sed():

    """
    This function ...
    :return:
    """

    # Determine the path to the SED file
    sed_path = fs.join(seds_path, "bcsed.dat")

    # Get the data
    wavelength_column, flux_column = np.loadtxt(sed_path, usecols=(0, 1), unpack=True)

    # Create the SED instance
    sed = SED()

    # Set the columns
    sed.table = tables.new([wavelength_column, flux_column], ["Wavelength", "Flux"])
    sed.table["Wavelength"].unit = "micron"
    sed.table["Flux"].unit = "W/m2" # = lambda * F_lambda !

    # Return the SED
    return sed

# -----------------------------------------------------------------

def load_example_zubko_sed():

    """
    This function ...
    :return:
    """

    # Determine the path to the SED file
    sed_path = fs.join(seds_path, "zubkosed.dat")

    # Get the data
    wavelength_column, flux_column = np.loadtxt(sed_path, usecols=(0, 2), unpack=True)

    # Create the SED instance
    sed = SED()

    # Set the columns
    sed.table = tables.new([wavelength_column, flux_column], ["Wavelength", "Flux"])
    sed.table["Wavelength"].unit = "micron"
    sed.table["Flux"].unit = "W/m2" # = lambda * F_lambda

    # Return the SED
    return sed

# -----------------------------------------------------------------

def load_example_themis_sed():

    """
    This function ...
    :return:
    """

    raise NotImplementedError("Not yet implemented")

    # Determine the path to the SED file
    # ...

# -----------------------------------------------------------------

max_stars_repo_name: wdobbels/CAAPR

max_issues_repo_name: wdobbels/CAAPR

max_forks_repo_name: wdobbels/CAAPR

reason: import numpy,from astropy


text:
import time
import os
from datetime import datetime

import numpy as np
from inflection import underscore

DASHES_TO_SPACES = {'_': ' ', '-': ' '}
REMOVE_PAREN = {'(': '', ')': ''}

def normalize_text(text, to_list=True, replace_chars = {'_': ' ', '-': ' ', '(': '', ')': ''}):
    text = underscore(text)  # converts to snake_case
    for old, new in replace_chars.items():
        text = text.replace(old, new)
    if to_list:
        return text.split(' ')
    else:
        return text

def unit_norm_rows(vectors):
    return vectors/np.linalg.norm(vectors, axis=1)[:,None]

def mean_of_rows(vectors):
    return np.mean(vectors, axis=0)

def max_of_rows(vectors):
    return np.max(vectors, axis=0)

def in_vocab(word_list, model):
    if isinstance(word_list, str):
        word_list = word_list.split(' ')
    return all([word in model.wv.vocab for word in word_list])

def get_dropped(old, new):
        return set(old).difference(set(new))

def no_op():
    return None

def path_to_name(path):
    return os.path.basename(path).split('.')[0]

def timeit(func, args=None):
    start = time.time()
    print('calling {0} \n'.format(func.__name__))
    result = func(*args) if args else func()
    print('{0} took {1} seconds \n\n'.format(func.__name__, time.time() - start))
    return result

def get_timestamp():
    return datetime.fromtimestamp(time.time()).strftime('%Y-%m-%d_%H-%M-%S')


text:
import numpy as np
from qsim import Operation, Circuit

class OperationZ(Operation):
    """
    Quantum X operation
    """
    def __init__(self, circuit: Circuit):
        op_matrix = np.array([
            [1, 0],
            [0, -1]])
        super().__init__(circuit, op_matrix)

max_stars_repo_name: ruanton/qsim

max_issues_repo_name: ruanton/qsim

max_forks_repo_name: ruanton/qsim


text:
"""
Reproduce Table 2

Condorcet Efficiencies under Random Society and
Spatial Model Assumptions
(201 voters, 5 candidates)

from

S. Merrill III, "A Comparison of Efficiency of Multicandidate
Electoral Systems", American Journal of Political Science, vol. 28,
no. 1, pp. 23-48, 1984.  :doi:`10.2307/2110786`

Typical result:

| Disp      |   1.0 |   1.0 |   1.0 |   1.0 |   0.5 |   0.5 |   0.5 |   0.5 |
| Corr      |   0.5 |   0.5 |   0.0 |   0.0 |   0.5 |   0.5 |   0.0 |   0.0 |
| Dims      |     2 |     4 |     2 |     4 |     2 |     4 |     2 |     4 |
|:----------|------:|------:|------:|------:|------:|------:|------:|------:|
| Plurality |  57.5 |  65.8 |  62.2 |  78.4 |  21.7 |  24.4 |  27.2 |  41.3 |
| Runoff    |  80.1 |  87.3 |  81.6 |  93.6 |  35.4 |  42.2 |  41.5 |  61.5 |
| Hare      |  79.2 |  86.7 |  84.0 |  95.4 |  35.9 |  46.8 |  41.0 |  69.9 |
| Approval  |  73.8 |  77.8 |  76.9 |  85.4 |  71.5 |  76.4 |  73.8 |  82.7 |
| Borda     |  87.1 |  89.3 |  88.2 |  92.3 |  83.7 |  86.3 |  85.2 |  89.4 |
| Coombs    |  97.8 |  97.3 |  97.9 |  98.2 |  93.5 |  92.3 |  93.8 |  94.5 |
| Black     | 100.0 | 100.0 | 100.0 | 100.0 | 100.0 | 100.0 | 100.0 | 100.0 |
| SU max    |  82.9 |  85.8 |  85.3 |  90.8 |  78.1 |  81.5 |  80.8 |  87.1 |
| CW        |  99.7 |  99.7 |  99.7 |  99.6 |  98.9 |  98.6 |  98.7 |  98.5 |

Many of these values match the paper closely, but some are consistently off by
up to 4%.
"""
import time
from collections import Counter
import numpy as np
from tabulate import tabulate
from elsim.methods import (fptp, runoff, irv, approval, borda, coombs,
                           black, utility_winner, condorcet)
from elsim.elections import normal_electorate, normed_dist_utilities
from elsim.strategies import honest_rankings, approval_optimal

n = 10_000
n_voters = 201
n_cands = 5

ranked_methods = {'Plurality': fptp, 'Runoff': runoff, 'Hare': irv,
                  'Borda': borda, 'Coombs': coombs, 'Black': black}

rated_methods = {'SU max': utility_winner,
                 'Approval': lambda utilities, tiebreaker:
                     approval(approval_optimal(utilities), tiebreaker)}

start_time = time.monotonic()

#             disp, corr, D
conditions = ((1.0, 0.5, 2),
              (1.0, 0.5, 4),
              (1.0, 0.0, 2),
              (1.0, 0.0, 4),
              (0.5, 0.5, 2),
              (0.5, 0.5, 4),
              (0.5, 0.0, 2),
              (0.5, 0.0, 4),
              )

results = []

for disp, corr, D in conditions:
    print(disp, corr, D)

    count = Counter()

    for iteration in range(n):
        v, c = normal_electorate(n_voters, n_cands, dims=D, corr=corr,
                                 disp=disp)

        """
        "Simulated utilities were normalized by range, that is, each voter's
        set of utilities were linearly expanded so that the highest and lowest
        utilities for each voter were 1 and 0, respectively."

        TODO: standard scores vs normalized don't matter for the ranked systems
        and don't affect approval much

        but This is necessary for the SU Maximizer results to match Merrill's.
        """
        utilities = normed_dist_utilities(v, c)
        rankings = honest_rankings(utilities)

        # If there is a Condorcet winner, analyze election, otherwise skip it
        CW = condorcet(rankings)
        if CW is not None:
            count['CW'] += 1

            for name, method in ranked_methods.items():
                if method(rankings, tiebreaker='random') == CW:
                    count[name] += 1

            for name, method in rated_methods.items():
                if method(utilities, tiebreaker='random') == CW:
                    count[name] += 1

    results.append(count)

elapsed_time = time.monotonic() - start_time
print('Elapsed:', time.strftime("%H:%M:%S", time.gmtime(elapsed_time)), '\n')

# Neither Tabulate nor Markdown support column span or multiple headers, but
# at least this prints to plain text in a readable way.
header = ['Disp\nCorr\nDims'] + [f'{x}\n{y}\n{z}' for x, y, z in conditions]

# Of those elections with CW, likelihood that method chooses CW
table = []
y_cw = np.array([c['CW'] for c in results])
for method in ('Plurality', 'Runoff', 'Hare', 'Approval', 'Borda', 'Coombs',
               'Black', 'SU max'):
    y = np.array([c[method] for c in results])
    table.append([method, *(y/y_cw*100)])

# Likelihood of Condorcet Winner (normalized by n iterations)
table.append(['CW', *(y_cw/n*100)])

print(tabulate(table, header, tablefmt="pipe", floatfmt='.1f'))

max_stars_repo_name: endolith/elsim

max_issues_repo_name: endolith/elsim

max_forks_repo_name: endolith/elsim


text:
#!/usr/bin/env python
# Wenchang Yang (wenchang@princeton.edu)
# Wed Jun 30 17:02:45 EDT 2021
if __name__ == '__main__':
    from misc.timer import Timer
    tt = Timer(f'start {__file__}')
import sys, os.path, os, glob, datetime
import xarray as xr, numpy as np, pandas as pd, matplotlib.pyplot as plt
#more imports
from scipy.optimize import minimize
from scipy.stats import genextreme
from numpy import exp, log
from tqdm import tqdm
from .shared import plot_smp_return_period, plot_return_period
#
if __name__ == '__main__':
    tt.check('end import')
#
#start from here
#negative log likelihood
def negLogLikelihood(params, data):
    """GEV negative log likelihood.
        params: (mu, sigma, xi)
        xx: sample(s)
    Ref: https://en.wikipedia.org/wiki/Generalized_extreme_value_distribution
    """
    mu, sigma, xi = params
    s = (data - mu)/sigma
    if xi == 0:
        return data.size*log(sigma) + np.sum(s + exp(-s))
    else:
        return -np.sum( log(xi*s>-1) ) + data.size*log(sigma) + np.sum( (1+1/xi)*log(1+xi*s) ) + np.sum( (1+xi*s)**(-1/xi) )
def fit(data, **kws):
    """GEV fit using maximum likelihood"""
    method = kws.pop('method', 'Nelder-Mead')
    xi_bounds = kws.pop('xi_bounds', (None, None))
    bounds = ( (None,None),  (0,None),  xi_bounds )
    x0_default = (data.mean(), data.std(), 0)
    x0 = kws.pop('x0', x0_default)

    return minimize(negLogLikelihood, x0, args=(data,), method=method, bounds=bounds, **kws)
def plot_fit(data, fit_kws=None, ax=None, **kws):
    if ax is None:
        fig,ax = plt.subplots()
    if fit_kws is None:
        fit_kws = {}
    upper = kws.pop('upper', data.size*10)
    #empirical return periods
    plot_smp_return_period(data, ax=ax, **kws)
    #fit by scipy.stats.genextreme
    c, mu, sigma = genextreme.fit(data)
    xi = -c
    print('genextreme fit:'.ljust(16), f'{mu=:.4g}; {sigma=:.4g};  {xi=:.4g}')
    #fit by user defined likelihood function
    r = fit(data, **fit_kws)
    if r.success:
        mu, sigma, xi = r.x
        print('wy fit:'.ljust(16), f'{mu=:.4g}; {sigma=:.4g}; {xi=:.4g}')
        plot_return_period(mu, sigma, xi, upper=upper, ax=ax, **kws)
    else:
        print(f'{r = }')
        print('[failed]:', r.message)
    ax.set_xlabel('return period')
    try:
        ax.set_ylabel(data.name)
    except:
        pass
    print()

    return r

def fit_bootstrap(data, nmc=100, mc_seed=None, **kws):
    """GEV fit bootstrap.
        data: input array-like data to fit
        nmc: size of Monte Carlo samples
        seed: np.random seed (default is 0)
    """
    rng = np.random.default_rng(mc_seed)
    datamc = np.random.choice(data, size=(nmc, data.size))
    params = np.zeros(shape=(nmc, 3)) + np.nan

    for ii in tqdm(range(nmc)):
        r = fit(datamc[ii,:], **kws)
        if r.success:
            params[ii,:] = r.x
        else:
            print(f'mc = {ii};', r.message)
    #with mp.Pool(processes=min(40, mp.cpu_count(), nmc)) as p:
    #    p.map(func_bs, range(nmc))
    mu = xr.DataArray(params[:,0], dims='mc')
    sigma = xr.DataArray(params[:,1], dims='mc')
    xi = xr.DataArray(params[:,2], dims='mc')

    ds = xr.Dataset(dict(mu=mu, sigma=sigma, xi=xi))
    return ds
def plot_fit_bootstrap(data, bsfit=None, nmc=100, mc_seed=None, ci=95, upper_rp=None, ax=None, fit_kws=None, **kws):
    if upper_rp is None:
        upper_rp = data.size*10
    if ax is None:
        fig,ax = plt.subplots()
    if fit_kws is None:
        fit_kws = {}
    #direct fit plot
    r = plot_fit(data, fit_kws=fit_kws, ax=ax, upper=upper_rp, **kws)
    mu, sigma, xi = r.x
    #bootstrap
    if bsfit is None: #do the bootstrap fit
        ds = fit_bootstrap(data, nmc=nmc, mc_seed=mc_seed, **fit_kws)
    else: #already done the bootstrap fit: use the result directly
        ds = bsfit
    ci_bounds = [(1-ci/100)/2, (1+ci/100)/2]
    for ii,daname in enumerate(('mu', 'sigma', 'xi')):
        q = ds[daname].quantile(ci_bounds, dim='mc')
        print(f'{daname} and {ci}% CI:'.rjust(20), f'{r.x[ii]:.4g}({q[0].item():.4g}, {q[1].item():.4g})')
    print()
    #confidence interval of the return value
    lower, upper = 1, np.log10(upper_rp) #return period bounds
    rp = np.logspace(lower, upper, 100)
    yy = [gev_return_period_inverse(rp, mu, sigma, xi)
        for mu,sigma,xi in zip(ds.mu.values, ds.sigma.values, ds.xi.values)]
    yy = xr.DataArray(yy, dims=('mc', 'rp')).assign_coords(rp=rp)
    yy.quantile(ci_bounds, dim='mc').plot(x='rp', ls='--', lw=1, hue='quantile', add_legend=False, **kws)

    ax.set_xlabel('return period')
    ax.set_ylabel('return value')

    return ds
def test(mu=None, sigma=None, xi=None, nmc=100, mc_seed=None, ci=95,nsmp=100):
    #specify params
    rng = np.random.default_rng()
    if mu is None:
        mu = rng.uniform(-10, 10)
    if sigma is None:
        sigma = rng.uniform(0, 10)
    if xi is None:
        xi = rng.uniform(-1, 1)
    true_values = mu,sigma,xi
    print('true params:'.ljust(16), f'{mu=:.4g}; {sigma=:.4g}; {xi=:.4g}')
    #generate data
    data = genextreme.rvs(-xi, loc=mu, scale=sigma, size=nsmp)
    #validate
    #fit_bootstrap plot
    fig,ax = plt.subplots()
    ds = plot_fit_bootstrap(data, nmc=nmc, mc_seed=mc_seed, ci=ci, ax=ax, color='C0')#, upper_rp=data.size*40)
    #fit summary
    r = fit(data)
    mu,sigma,xi = r.x
    ci_bounds = [(1-ci/100)/2, (1+ci/100)/2]
    s = ''
    danames = ('mu', 'sigma', 'xi')
    pnames = ('$\\mu$:', '$\\sigma$:', '$\\xi$:')
    for ii,(daname,pname, tv) in enumerate(zip(danames, pnames, true_values)):
        q = ds[daname].quantile(ci_bounds, dim='mc')
        s += pname + f' {r.x[ii]:.4g}(true:{tv:.4g}); {ci}% CI: ({q[0].item():.4g}, {q[1].item():.4g})\n'
    ax.text(1, 0, s, transform=ax.transAxes, ha='right', va='bottom', fontsize='small', alpha=0.5)
    print(s)

if __name__ == '__main__':
    from wyconfig import * #my plot settings
    plt.close('all')
    if len(sys.argv)<=1:
        test()
    elif len(sys.argv)>1 and sys.argv[1]=='test': #e.g. python -m wyextreme.gev test xi=-0.1
        kws = dict(mu=None, sigma=None, xi=None, nmc=100, mc_seed=None, ci=95, nsmp=100)
        if len(sys.argv)>2:
            for s in sys.argv[2:]:
                key,v = s.split('=')
                v = int(v) if key in ('nmc', 'mc_seed', 'nsmp') else float(v)
                if key in kws: kws[key] = v
        test(**kws)
    elif len(sys.argv)>2: # two input data files to compare
        da0 = xr.open_dataarray(sys.argv[1])
        da1 = xr.open_dataarray(sys.argv[2])
        fig, ax = plt.subplots()
        plot_fit_bootstrap(da0, ax=ax, color='C0')
        plot_fit_bootstrap(da1, ax=ax, color='C1')
    elif len(sys.argv)>1:#read input data file and do the analysis; otherwise
        ifile = sys.argv[1]
        da = xr.open_dataarray(ifile)
        plot_fit(da)

    #savefig
    if 'savefig' in sys.argv:
        figname = __file__.replace('.py', f'.png')
        wysavefig(figname)
    tt.check(f'**Done**')
    plt.show()

reason: from numpy,from scipy


text:
import numpy as np

# ============================================================
# Basic Functions
# ============================================================

def _assign_interval_base(x, boundaries):
    """Assign each value in x an interval from boundaries.

    Parameters
    ----------
    x: numpy.array, shape (number of examples,)
        The column of data that need to be discretized.

    boundaries: numpy.array, shape (number of interval boundaries,)
        The boundary values of the intervals to discretize target x.

    Returns
    -------
    intervals: numpy.ndarray, shape (number of examples,2)
        The array of intervals that are closed to the right.
        The left column and right column of the array are the
        left and right boundary respectively.
    """
    # Add -inf and inf to the start and end of boundaries
    max_value = max(x)
    boundaries = np.unique(
        np.concatenate((np.array([-float('inf')]),
                        boundaries,np.array([float('inf')])),
                        axis=0))
    # The max boundary that is smaller than x_i is its lower boundary.
    # The min boundary that is >= than x_i is its upper boundary.
    # Adding equal is because all intervals here are closed to the right.
    boundaries_diff_boolean = x.reshape(1,-1).T > boundaries.reshape(1,-1)
    lowers = np.array([boundaries[b].max() for b in boundaries_diff_boolean])
    uppers = np.array([boundaries[b].min() for b in ~boundaries_diff_boolean])

    # Replace the upper value with inf if the value equals the maximum feature value
    n = x.shape[0]
    uppers = np.where(uppers==max_value,[float('inf')]*n,uppers)
    # Array of intervals that are closed to the right
    intervals= np.stack((lowers, uppers), axis=1)
    return intervals

def assign_interval_unique(x, boundaries):
    """Assign each value in x an interval from boundaries.

    Parameters
    ----------
    x: numpy.array, shape (number of examples,)
        The column of data that need to be discretized.

    boundaries: numpy.array, shape (number of interval boundaries,)
        The boundary values of the intervals to discretize target x.

    Returns
    -------
    intervals: numpy.ndarray, shape (number of examples,2)
        The array of intervals that are closed to the right.
        The left column and right column of the array are the
        left and right boundary respectively.

    unique_intervals: numpy.ndarray, shape (number of unique intervals,2)
        The unique intervals that are closed to the right.
        The left column and right column of the array are the
        left and right boundary respectively.
    """
    intervals= _assign_interval_base(x, boundaries)
    unique_intervals = np.unique(intervals, axis=0)
    return intervals, unique_intervals

def assign_interval_str(x, boundaries, delimiter='~'):
    """Assign each value in x an interval from boundaries.

    Parameters
    ----------
    x: numpy.array, shape (number of examples,)
        The column of data that need to be discretized.

    boundaries: numpy.array, shape (number of interval boundaries,)
        The boundary values of the intervals to discretize target x.

    delimiter: string, optional(default='~')
        The returned array will be an array of intervals. Each interval is
        representated by string (i.e. '1~2'), which takes the form
        lower+delimiter+upper. This parameter control the symbol that
        connects the lower and upper boundaries.
    Returns
    -------
    intervals_str: numpy.array, shape (number of examples,)
        Discretized result. The array of intervals that represented by strings.
    """
    intervals= _assign_interval_base(x, boundaries)
    # use join rather than use a+delimiter+b makes this line faster
    intervals_str = np.array(
        [delimiter.join((str(a),str(b))) for a,b in zip(intervals[:,0],
                                                        intervals[:,1])]
        )
    return intervals_str

def pivot_table_np(index, column):
    """Perform cross-tabulation to index vector and column vector.
    The returned pivot table has unique index values as its rows,
    unique column values as its columns, and sample frequencies
    of index-column combinations as its values.

    Parameters
    ----------
    index: numpy.array, shape (number of examples,)
             The vector whose unique values will be the rows of pivot table

    column: numpy.array, shape (number of examples,)
            The vector whose unique values will be the columns of pivot table

    Returns
    -------
    pivot_table: numpy.array, shape (number of unqiue index values,
                                     number of unqiue column values)
             The pivot table that has unique index values as its rows,
             unique column values as its columns, and sample frequencies
             of index-column combinations as its values.

    column_unique: numpy.array, shape (number of unqiue column values,)
            Column labels of the pivot table.
            The unqiue column values are sorted ascendingly.

    index_unique: numpy.array, shape (number of unqiue index values,)
            Row labels of the pivot table.
            The unqiue index values are sorted ascendingly.
    """

    # unique values of index and column
    column_unique, index_unique = np.unique(column), np.unique(index)
    n_j = column_unique.shape[0] # number of columns, number of indexes (rows)

    # Sample frequencies of each combination of index vector and column vector
    # Consider the value array has two hierarchical indexes (index & column)
    dot_multiply = np.stack((index, column), axis=1)
    dot_multiply, value = np.unique(dot_multiply, axis=0, return_counts=True)

    # generate pivot table, there may be duplicated rows for index vector
    zeros = np.zeros(dot_multiply.shape[0])
    pivot_table = np.vstack(
        [np.where(dot_multiply[:,1]==j, value, zeros) for j in column_unique]
        ).T

    # merge the pivot table rows with the same index value
    lowers = dot_multiply[:,0]
    pivot_table = np.vstack(
        [pivot_table[lowers==i,:].sum(axis=0) for i in index_unique]
        )
    return pivot_table, column_unique, index_unique

def interval_to_boundary_vector(vector, delimiter='~'):
    """Transform an array of interval strings into the
    unique boundaries of such intervals.

    Parameters
    ----------
    vector: numpy.array, shape (number of examples,)
        The array of interval whose unique boundaries will
        be returned.

    delimiter: string, optional(default='~')
        The interval is representated by string (i.e. '1~2'),
        which takes the form lower+delimiter+upper. This parameter
        control the symbol that connects the lower and upper boundaries.
    Returns
    -------
    boundaries: numpy.array, shape (number of interval boundaries,)
        An array of boundary values.
    """
    boundaries = np.array(list(set(delimiter.join(np.unique(vector)).split(delimiter))))
    boundaries = boundaries[(boundaries!='-inf') & (boundaries!='inf')].astype(float)
    return boundaries

def map_np(array, dictionary):
    """map function for numpy array
    Parameters
    ----------
    array: numpy.array, shape (number of examples,)
            The array of data to map values to.

    distionary: dict
            The distionary object.

    Return
    ----------
    result: numpy.array, shape (number of examples,)
            The mapped result.
    """
    return [dictionary[e] for e in array]


text:
import tensorflow as tf
import numpy as np

# DISCLAIMER:
# Parts of this code file are derived from
# https://github.com/williamleif/GraphSAGE
# which is under an identical MIT license as CAGNIR

def uniform(shape, scale=0.05, dtype=tf.float32, name=None):
    """Uniform init."""
    initial = tf.random_uniform(shape, minval=-scale, maxval=scale, dtype=dtype)
    return tf.Variable(initial, name=name)

def glorot(shape, get_var=False, name=None, dtype=tf.float32):
    """Glorot & Bengio (AISTATS 2010) init."""

    init_range = np.sqrt(6.0/(int(shape[0])+shape[1]))

    if get_var:

        if name is None:
            raise Exception('name is not defined.')

        initial = tf.contrib.layers.xavier_initializer()
        return tf.get_variable(name=name, shape=shape, initializer=initial, dtype=dtype)

    else:

        initial = tf.contrib.layers.xavier_initializer()
        return tf.Variable(initial, dtype=dtype, name=name)

def zeros(shape, get_var=False, name=None, dtype=tf.float32):
    """All zeros."""
    if get_var:

        if name is None:
            raise Exception('name is not defined.')

        return tf.get_variable(name=name, shape=shape, initializer=tf.zeros_initializer(), dtype=dtype)

    else:

        initial = tf.zeros(shape, dtype=dtype)
        return tf.Variable(initial, name=name)

def ones(shape, name=None, dtype=tf.float32):
    """All ones."""
    initial = tf.ones(shape, dtype=dtype)
    return tf.Variable(initial, name=name)

max_stars_repo_name: rmhsiao/CAGNIR

max_issues_repo_name: rmhsiao/CAGNIR

max_forks_repo_name: rmhsiao/CAGNIR


text:
import numpy as np
import cv2

class DarknetUtil:
    DARKNET_SCALE = 0.00392
    def __init__(self, cfg_path : str, model_path : str, net_size : tuple = (416, 416)):
        self.net = cv2.dnn.readNetFromDarknet(cfg_path, model_path)
        self.width = net_size[0]
        self.height = net_size[1]

    def __del__(self):
        pass

    def detect(self, image, thresh = 0.25):
        Height, Width = image.shape[:2]
        blob = cv2.dnn.blobFromImage(image, self.DARKNET_SCALE, (self.width, self.height), (0,0,0), True, crop=False)

        self.net.setInput(blob)
        layer_names = self.net.getLayerNames()
        output_layers = [layer_names[i - 1] for i in self.net.getUnconnectedOutLayers()]
        outs = self.net.forward(output_layers)

        class_ids = []
        confidences = []
        boxes = []
        conf_threshold = 0.5
        nms_threshold = 0.4

        for out in outs:
            for detection in out:
                scores = detection[5:]
                class_id = np.argmax(scores)
                confidence = scores[class_id]
                if confidence > thresh:
                    center_x = int(detection[0] * Width)
                    center_y = int(detection[1] * Height)
                    w = int(detection[2] * Width)
                    h = int(detection[3] * Height)
                    x = int(center_x - w / 2)
                    y = int(center_y - h / 2)
                    class_ids.append(class_id)
                    confidences.append(float(confidence))
                    boxes.append([x, y, w, h])

        indices = cv2.dnn.NMSBoxes(boxes, confidences, conf_threshold, nms_threshold)

        retBoxes = []
        for i in indices:
            x,y,w,h = boxes[i]
            retBoxes.append({"x":x, "y":y, "w":w, "h":h, "id":class_ids[i], "thresh":confidences[i]})

        return retBoxes

if __name__ == "__main__":
    pass


text:
import random
import numpy as np
import torch
from torch import Tensor
from typing import Tuple
import torch.nn as nn
from transformers import AutoModel
import torch.nn.functional as F

"""
- DATA: we split each docs to 2 parts: 75% and 25%
- Model: we formulate this problem as seq2seq and mask prediction task, model consist of 2 components: encoder and decoder
    + Encoder, we use 2 bi-gru for learning local and global features
    + Decoder has 2 parts: the first reconstruct invert 75% docs (.i.e from token T to 0), the rest will be mask prediction with 50%
"""

class Encoder(nn.Module):
    def __init__(self,
                 hidden_size,
                 device,
                 weight_top_embedding,
                 weight_bot_embedding,
                 num_layers=2,
                 bidirectional=True,
                 p=0.2,
                 batch_first=True):

        super(Encoder, self).__init__()

        self.num_layers = num_layers
        self.hidden_size = hidden_size
        self.D = 2 if bidirectional else 1
        self.bidirectional = bidirectional
        self.device = device
        self.top_embedd_layer = nn.Embedding.from_pretrained(embeddings=weight_top_embedding.vectors)
        self.bot_embedd_layer = nn.Embedding.from_pretrained(embeddings=weight_bot_embedding.vectors)

        # the first rnn to learning sentiment text, using sentiment word embedding
        self.rnn1 = nn.GRU(input_size=weight_top_embedding.dim, hidden_size=hidden_size, num_layers=num_layers,
                           dropout=p, bidirectional=bidirectional, batch_first=batch_first)

        # the second rnn to learning global represent of text, using fasttex
        self.rnn2 = nn.GRU(input_size=weight_bot_embedding.dim, hidden_size=hidden_size, num_layers=num_layers,
                           dropout=p, bidirectional=bidirectional, batch_first=batch_first)

        self.dropout = nn.Dropout(p)

    def _init_hidden(self, batch_size):
        """
        Args:
            batch_size: int
        Returns:
            h (Tensor): the first hidden state of rnn
        """

        if self.bidirectional:
            h = torch.zeros((self.num_layers*2, batch_size, self.hidden_size))
        else:
            h = torch.zeros((self.num_layers, batch_size, self.hidden_size))
        nn.init.xavier_normal_(h)

        return h.to(self.device)

    def forward(self, top_x, bot_x):
        """
        Args:
            top_x, bot_x: Tensor(batch_size, seq_len)
        Returns:
            ouput1, h1: Tuple(Tensor, Tensor)
            output2, h2: Tuple(Tensor, Tensor)
        """

        top_embedded = self.top_embedd_layer(top_x)
        bot_embedded = self.bot_embedd_layer(bot_x)

        h = self._init_hidden(top_x.size(0))

        output1, h1 = self.rnn1(top_embedded, h) # h1, h2 size of (D ∗ num_layers, batch_size, hidden_size)
        output2, h2 = self.rnn2(bot_embedded, h) # output1, output2 size of (batch_size, seq_len, D * hidden_size)

        return (output1, h1), (output2, h2)

class BahdanauAttention(nn.Module):
    def __init__(self, dec_dim: int, enc_dim: int, num_hiddens: int):
        super().__init__()
        self.W1 = nn.Linear(dec_dim, num_hiddens, bias=False)
        self.W2 = nn.Linear(enc_dim, num_hiddens, bias=False)
        self.v = nn.Linear(num_hiddens, 1, False)

    def forward(self, query: Tensor, value: Tensor) -> Tuple[Tensor, Tensor]:
        """
        Args:
            value (Tensor(batch_size, seq_len, encoder hidden dimension): the hidden_state of tokens in encoder
            query (Tensor(batch_size, 1, decoder hidden dimension)): the hidden state of decoder at time step t
        Returns:
            attention_weight (Tensor)
            context_vector (Tensor)
        """

        score = self.v(torch.tanh(self.W1(query) + self.W2(value))) # `score` size of: (batch_size, seq_len, 1)

        attention_weight = F.softmax(score.squeeze(-1), dim=1) # `attention` size of: (batch_size, seq_len)

        context_vector = torch.bmm(attention_weight.unsqueeze(1), value).squeeze(1) # `context_vector` size of: (batch_size, num_hiddens)
        return attention_weight, context_vector

class Decoder(nn.Module):
    def __init__(self,
                 hidden_size,
                 weight_top_embedding,
                 weight_bot_embedding,
                 encoder_output_dim,
                 num_layers=2,
                 bidirectional=True,
                 p=0.2,
                 batch_first=True):

        super(Decoder, self).__init__()

        self.num_layers = num_layers
        self.hidden_size = hidden_size
        self.D = 2 if bidirectional else 1
        self.bidirectional = bidirectional
        self.vocab_size = weight_top_embedding.vectors.size(0)

        self.top_embedd_layer = nn.Embedding.from_pretrained(embeddings=weight_top_embedding.vectors)
        self.bot_embedd_layer = nn.Embedding.from_pretrained(embeddings=weight_bot_embedding.vectors)

        self.rnn1 = nn.GRU(input_size=weight_top_embedding.dim, hidden_size=hidden_size, num_layers=num_layers,
                           dropout=p, bidirectional=bidirectional, batch_first=batch_first)

        self.rnn2 = nn.GRU(input_size=weight_bot_embedding.dim, hidden_size=hidden_size, num_layers=num_layers,
                           dropout=p, bidirectional=bidirectional, batch_first=batch_first)

        self.dropout = nn.Dropout(p)

        self.attention1 = BahdanauAttention(dec_dim=hidden_size * self.D, enc_dim=encoder_output_dim,
                                            num_hiddens=hidden_size * self.D)
        self.attention2 = BahdanauAttention(dec_dim=hidden_size * self.D, enc_dim=encoder_output_dim,
                                            num_hiddens=hidden_size * self.D)

        self.fc = nn.Sequential(nn.Linear(2*(encoder_output_dim + hidden_size*2), 2 * encoder_output_dim),
                                nn.ReLU())

    def forward(self, x1: Tensor, output1_encoder: Tensor, h1_encoder: Tensor,
                x2: Tensor, output2_encoder: Tensor, h2_encoder: Tensor) -> Tuple[Tensor, Tensor, Tensor]:

        """
        Args:
            `x1`, `x2` (Tensor): seq number for embedding layer
            `output1_encoder`, `h1_encoder`, `output2_encoder`, `h2_encoder` (Tensor): output of encoder layer
        """

        # `x1`, `x2` size of batch_size -> (batch_size, 1): `1` denote seq_length
        x1, x2 = x1.unsqueeze(1), x2.unsqueeze(1)

        top_embedd = self.dropout(self.top_embedd_layer(x1))
        bot_embedd = self.dropout(self.bot_embedd_layer(x2)) # (batch_size, 1, embedd dim)

        output1, h1 = self.rnn1(top_embedd, h1_encoder)
        output2, h2 = self.rnn2(bot_embedd, h2_encoder)

        _, context_vector1 = self.attention1(output1, output1_encoder)
        _, context_vector2 = self.attention2(output2, output2_encoder)

        output1 = torch.cat((output1.squeeze(1), context_vector1), dim=1)
        output2 = torch.cat((output2.squeeze(1), context_vector2), dim=1)
        output = torch.cat((output1, output2), dim=1) # [batch_size, 2* (encoder_output_dim + hidden_size*2)]
        output = self.fc(output) # [batch_size, 2 * encoder_output_dim]

        return output, h1, h2

class SentimentModel(nn.Module):
    def __init__(self, encoder, decoder1, decoder2, num_aspect):
        super(SentimentModel, self).__init__()
        self.encoder = encoder
        self.decoder1 = decoder1
        self.decoder2 = decoder2
        self.relu = nn.ReLU()
        self.T = nn.Linear(2*(encoder.hidden_size*encoder.D), num_aspect)
        self.Z1 = nn.Linear(num_aspect, decoder1.vocab_size)
        self.Z2 = nn.Linear(num_aspect, decoder2.vocab_size)

    def forward(self, top_source, bot_source, top_res, top_mask, teacher_forcing_ratio=0.5):
        """
        Args:
            `top_soucre`, `bot_source` (Tensor): input of encoder block
            `top_res`, `bot_res` (Tensor): input of the first decoder block
            `top_mask`, `bot_mask` (Tensor): input of the second decoder block
        """
        batch_size = top_mask.size(0)
        mask_len = top_mask.size(1)
        res_len = top_res.size(1)

        output_reverts = torch.zeros(res_len, batch_size, self.decoder1.vocab_size)
        output_masks = torch.zeros(mask_len, batch_size, self.decoder2.vocab_size)

        (output1, h1), (output2, h2) = self.encoder(top_source, bot_source)

        input_ = top_res[:, 0]
        for i in range(1, res_len):
            output_revert, h1, h2 = self.decoder1(x1=input_, output1_encoder=output1, h1_encoder=h1, x2=input_, output2_encoder=output2, h2_encoder=h2)

            output_revert = self.relu(self.T(output_revert))
            output_revert = self.Z1(output_revert)

            output_reverts[i] = output_revert
            teacher_force = random.random() < teacher_forcing_ratio

            top1 = output_revert.argmax(1)
            input_ = top_res[:, i] if teacher_force else top1

        (output1, h1), (output2, h2) = self.encoder(top_source, bot_source)

        input_ = top_mask[:, 0]
        for i in range(1, mask_len):
            output_mask, h1, h2 = self.decoder2(x1=input_, output1_encoder=output1, h1_encoder=h1, x2=input_, output2_encoder=output2, h2_encoder=h2)

            output_mask = self.relu(self.T(output_mask))
            output_mask = self.Z2(output_mask)

            output_masks[i] = output_mask
            teacher_force1 = random.random() < teacher_forcing_ratio

            top1 = output_mask.argmax(1)
            input_ = top_mask[:, i] if teacher_force1 else top1

        # `output_reverts`, `output_mask` size of : (revert seq_len/mask seq_len, batch_size, decoder1.vocab_size/decoder2.vocab_size)
        return output_reverts, output_masks

class EarlyStopping:
    def __init__(self, patience=5, delta=0):
        self.patience = patience
        self.counter = 0
        self.best_score = None
        self.early_stop = False
        self.val_loss_min = np.Inf
        self.delta = delta

    def __call__(self, val_loss):
        score = -val_loss
        if self.best_score is None:
            self.best_score = score
        elif score < self.best_score + self.delta:
            self.counter += 1
            if self.counter >= self.patience:
                self.early_stop = True
        else:
            self.best_score = score
            self.counter = 0

max_stars_repo_name: manhtrantienhn/sentiment

max_issues_repo_name: manhtrantienhn/sentiment

max_forks_repo_name: manhtrantienhn/sentiment


text:
# _*_ codig utf8 _*_

import numpy as np
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import matplotlib as mpl

mpl.rcParams['font.sans-serif'] = [u'simHei']

def f(x,y):# 原函数
    return x**2 + y ** 2

def h(x):# 导函数
    return 2*x

X=[]
Y=[]
Z=[]
x=2
y=2
step=0.1
f_change=f(x,y)
f_current=f(x,y)
X.append(x)
Y.append(y)
Z.append(f_current)
while f_change > 1e-10:
    x = x - step * h(x)
    y=y-step*h(y)
    f_change = f_current-f(x,y)
    f_current = f(x,y)
    X.append(x)
    Y.append(y)
    Z.append(f_current)

print('结果为(%s,%s)' %(x,y))
print(X)
fig = plt.figure()
ax=Axes3D(fig)
X2 = np.arange(-2,2,0.2)
Y2 = np.arange(-2,2,0.2)
X2,Y2=np.meshgrid(X2,Y2)
Z2=X2**2+Y2**2
ax.plot_surface(X2,Y2,Z2,rstride=1,cstride=1,cmap='rainbow')
ax.plot(X,Y,Z,'bo--')
#plt.plot(X2,Y2,'-',color='#666666')
#plt.plot(X,Y,'bo--')
ax.set_title('$y=x^2+y^2$函数求解最小值，最终解为（%.2f,%.2f,%.2f）'%(x,y,f_current))
plt.show()


text:
import numpy as np
import time

def euclidian_distance(x, y):
    d = x-y
    return np.sqrt(np.dot(d, d))

def manhattan_distance(x, y):
    return np.sum(np.abs(x-y))

#@profile
def fill_dist_array(x, y, d_f):
    tmp_d_a = np.ndarray((x.shape[0], y.shape[0]))
    for i in xrange(x.shape[0]):
        for j in xrange(y.shape[0]):
            tmp_d_a[i,j] = d_f(i,j)
    return tmp_d_a

#@profile
def DTW(x, y, dist_function=None, dist_array=None):
    """
    Default is euclidian distance, otherwise provide (in order of priority):
     - a distance function as dist_function
     - a distance array as dist_array[x_ind][y_ind]
    """

    d_f = None
    d_a = dist_array
    if dist_function == None and dist_array == None:
        d_f = lambda i,j: euclidian_distance(x[i], y[j])
    else:
        if dist_function != None:
            d_f = lambda i,j: dist_function(x[i], y[j])
        elif dist_array != None:
            d_f = lambda i,j: dist_array[i][j]
    #if d_a == None:
    #    d_a = fill_dist_array(x, y, d_f)

    cost = np.empty((x.shape[0], y.shape[0]), dtype=np.float)
    cost[0,0] = d_f(0,0)
    N = x.shape[0]
    for i in xrange(1, N):
        cost[i,0] = d_f(i,0)
    M = y.shape[0]
    for j in xrange(1, M):
        cost[0,j] = d_f(0,j)
    # the dynamic programming loop
    for i in xrange(1, N):
        for j in xrange(1, M):
            #cost[i,j] = d_a[i,j] + np.min([cost[i-1,j], cost[i-1,j-1],
            cost[i,j] = d_f(i,j) + np.min((cost[i-1,j], cost[i-1,j-1],
                                                         cost[i,j-1]))
    # here, cost[x.shape[0]-1,y.shape[0]-1] is the best possible distance
    path_x = None
    path_y = None
    return cost[N-1,M-1], cost, (path_x, path_y)

def test():
    np.random.seed(42)
    a = np.random.random((170, 10))
    b = np.random.random((130, 10))
    t = time.time()
    for k in xrange(2):
        d = DTW(a, b, euclidian_distance)
    print d
    print "took:", ((time.time() - t) / k),  "seconds per run"

if __name__ == '__main__':
    test()


text:
import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
from sklearn.datasets import load_digits
from scipy.interpolate import interp1d
from scipy.interpolate import interp2d
from sklearn.decomposition import PCA
from sklearn import datasets
import matplotlib.cm as cm
import cv2
import matplotlib.image as mpimg
from sklearn import svm
from sklearn import decomposition
from pdb import set_trace

def interpol_im(im, dim1 = 8, dim2 = 8, plot_new_im = False, cmap = 'binary', grid_off = False):

	if (len(im.shape) == 3):
		im = im[:,:,0]

	x = np.arange(im.shape[1])
	y = np.arange(im.shape[0])

	f2d = interp2d(x,y,im)

	x_new = np.linspace(0, im.shape[1], dim1)
	y_new = np.linspace(0, im.shape[0], dim2)

	im_new = f2d(x_new,y_new)

	if (plot_new_im):

		plt.imshow(im_new,interpolation='nearest',cmap=cmap)
		if (grid_off):
			plt.grid('off')
		plt.show()

	return im_new.flatten()

def pca_svm_pred(imfile, md_pca, md_clf, dim1 = 45, dim2 = 60):

	flatten_im = interpol_im(imfile, dim1 = dim1, dim2 = dim2)

	flatten_im = flatten_im.reshape(1,-1)

	flatten_im_proj = md_pca.transform(flatten_im)

	predict = md_clf.predict(flatten_im_proj)

	return predict[0]

def pca_X(X, n_comp = 10):

	md_pca = PCA(n_comp,whiten = True)
	md_pca.fit(X)
	X_proj = md_pca.transform(X)

	return md_pca, X_proj

def rescale_pixel(unseen):

	im_flat = interpol_im(unseen)

	im_flat = im_flat * 15
	im = im_flat.astype(int)

	im = 15 - im

	return im

def svm_train(X, y, gamma = 0.001, C = 100):

	md_clf = svm.SVC(gamma = gamma, C = C)
	md_clf.fit(X, y)

	return md_clf

reason: import numpy,from scipy


text:
# -*- coding: utf-8 -*-
#
# Authors: Alexandre Gramfort <alexandre.gramfort@inria.fr>
#          Eric Larson <larson.eric.d@gmail.com>
#          Joan Massich <mailsik@gmail.com>
#          Guillaume Favelier <guillaume.favelier@gmail.com>
#          Oleh Kozynets <ok7mailbox@gmail.com>
#
# License: Simplified BSD

import pytest
import numpy as np
import os.path as path
from mne import read_source_estimate
from mne.datasets import testing
from mne.viz._brain import _Brain
from mne.viz._brain.colormap import _calculate_lut

from matplotlib import cm

data_path = testing.data_path(download=False)
subject_id = 'sample'
subjects_dir = path.join(data_path, 'subjects')
fname_stc = path.join(data_path, 'MEG/sample/sample_audvis_trunc-meg')
surf = 'inflated'

@testing.requires_testing_data
def test_brain_init(renderer):
    """Test initialization of the _Brain instance."""
    backend_name = renderer.get_3d_backend()
    hemi = 'both'

    with pytest.raises(ValueError, match='hemi'):
        _Brain(subject_id=subject_id, hemi="split", surf=surf)
    with pytest.raises(ValueError, match='figure'):
        _Brain(subject_id=subject_id, hemi=hemi, surf=surf, figure=0)
    with pytest.raises(ValueError, match='interaction'):
        _Brain(subject_id=subject_id, hemi=hemi, surf=surf, interaction=0)
    with pytest.raises(KeyError):
        _Brain(subject_id=subject_id, hemi="foo", surf=surf)

    brain = _Brain(subject_id, hemi, surf, subjects_dir=subjects_dir)
    if backend_name != 'mayavi':
        brain.show()

@testing.requires_testing_data
def test_brain_add_data(renderer):
    """Test adding data in _Brain instance."""
    backend_name = renderer.get_3d_backend()
    stc = read_source_estimate(fname_stc)

    hemi = 'lh'
    hemi_data = stc.data[:len(stc.vertices[0]), 10]
    hemi_vertices = stc.vertices[0]
    fmin = stc.data.min()
    fmax = stc.data.max()

    brain_data = _Brain(subject_id, hemi, surf, size=300,
                        subjects_dir=subjects_dir)

    with pytest.raises(ValueError):
        brain_data.add_data(array=np.array([0, 1, 2]))
    with pytest.raises(ValueError):
        brain_data.add_data(hemi_data, fmin=fmin, hemi=hemi,
                            fmax=fmax, vertices=None)

    brain_data.add_data(hemi_data, fmin=fmin, hemi=hemi, fmax=fmax,
                        colormap='hot', vertices=hemi_vertices,
                        colorbar=False)

    if backend_name != 'mayavi':
        brain_data.show()

def test_brain_colormap():
    """Test brain's colormap functions."""
    colormap = "coolwarm"
    alpha = 1.0
    fmin = 0.0
    fmid = 0.5
    fmax = 1.0
    center = None
    _calculate_lut(colormap, alpha=alpha, fmin=fmin,
                   fmid=fmid, fmax=fmax, center=center)
    center = 0.0
    colormap = cm.get_cmap(colormap)
    _calculate_lut(colormap, alpha=alpha, fmin=fmin,
                   fmid=fmid, fmax=fmax, center=center)

max_stars_repo_name: fmamashli/mne-python

max_issues_repo_name: fmamashli/mne-python

max_forks_repo_name: fmamashli/mne-python


text:
# -*- coding: utf-8 -*-
import os
import sys
import argparse
import datetime
import pdb
import logging
import pprint
import numpy
import torch

from eznlp import auto_device
from eznlp.token import LexiconTokenizer
from eznlp.dataset import Dataset
from eznlp.config import ConfigDict
from eznlp.model import OneHotConfig, MultiHotConfig, EncoderConfig, CharConfig, SoftLexiconConfig
from eznlp.model import ELMoConfig, BertLikeConfig, FlairConfig
from eznlp.model import SequenceTaggingDecoderConfig, SpanClassificationDecoderConfig, BoundarySelectionDecoderConfig
from eznlp.model import ExtractorConfig
from eznlp.model.bert_like import segment_uniformly_for_bert_like
from eznlp.training import Trainer, count_params, evaluate_entity_recognition

from utils import add_base_arguments, parse_to_args
from utils import load_data, dataset2language, load_pretrained, load_vectors, build_trainer, header_format, profile

def parse_arguments(parser: argparse.ArgumentParser):
    parser = add_base_arguments(parser)

    group_data = parser.add_argument_group('dataset')
    group_data.add_argument('--dataset', type=str, default='conll2003',
                            help="dataset name")
    group_data.add_argument('--doc_level', default=False, action='store_true',
                            help="whether to load data at document level")
    group_data.add_argument('--corrupt_rate', type=float, default=0.0,
                            help="boundary corrupt rate")
    group_data.add_argument('--pipeline', default=False, action='store_true',
                            help="whether to save predicted chunks for pipeline")

    group_decoder = parser.add_argument_group('decoder configurations')
    group_decoder.add_argument('--ck_decoder', type=str, default='sequence_tagging',
                               help="chunk decoding method", choices=['sequence_tagging', 'span_classification', 'boundary_selection'])
    # Loss
    group_decoder.add_argument('--fl_gamma', type=float, default=0.0,
                               help="Focal Loss gamma")
    group_decoder.add_argument('--sl_epsilon', type=float, default=0.0,
                               help="Label smoothing loss epsilon")

    # Sequence tagging
    group_decoder.add_argument('--scheme', type=str, default='BIOES',
                               help="sequence tagging scheme", choices=['BIOES', 'BIO2'])
    group_decoder.add_argument('--no_crf', dest='use_crf', default=True, action='store_false',
                               help="whether to use CRF")

    # Span-based
    group_decoder.add_argument('--agg_mode', type=str, default='max_pooling',
                               help="aggregating mode")
    group_decoder.add_argument('--num_neg_chunks', type=int, default=100,
                               help="number of sampling negative chunks")
    group_decoder.add_argument('--max_span_size', type=int, default=10,
                               help="maximum span size")
    group_decoder.add_argument('--ck_size_emb_dim', type=int, default=25,
                               help="span size embedding dim")

    # Boundary selection
    group_decoder.add_argument('--no_biaffine', dest='use_biaffine', default=True, action='store_false',
                               help="whether to use biaffine")
    group_decoder.add_argument('--affine_arch', type=str, default='FFN',
                               help="affine encoder architecture")
    group_decoder.add_argument('--neg_sampling_rate', type=float, default=1.0,
                               help="Negative sampling rate")
    group_decoder.add_argument('--hard_neg_sampling_rate', type=float, default=1.0,
                               help="Hard negative sampling rate")
    group_decoder.add_argument('--hard_neg_sampling_size', type=int, default=5,
                               help="Hard negative sampling window size")
    group_decoder.add_argument('--sb_epsilon', type=float, default=0.0,
                               help="Boundary smoothing loss epsilon")
    group_decoder.add_argument('--sb_size', type=int, default=1,
                               help="Boundary smoothing window size")
    group_decoder.add_argument('--sb_adj_factor', type=float, default=1.0,
                               help="Boundary smoothing probability adjust factor")
    return parse_to_args(parser)

def collect_IE_assembly_config(args: argparse.Namespace):
    drop_rates = (0.0, 0.05, args.drop_rate) if args.use_locked_drop else (args.drop_rate, 0.0, 0.0)

    ohots_config = {}
    if args.emb_dim > 0:
        vectors = load_vectors(args.language, args.emb_dim, unigram=True)
        ohots_config['text'] = OneHotConfig(field='text', vectors=vectors, emb_dim=args.emb_dim, freeze=args.emb_freeze)

    if args.language.lower() == 'chinese' and args.use_bigram:
        vectors = load_vectors(args.language, 50, bigram=True)
        ohots_config['bigram'] = OneHotConfig(field='bigram', vectors=vectors, emb_dim=50, freeze=args.emb_freeze)

    ohots_config = ConfigDict(ohots_config) if len(ohots_config) > 0 else None

    if args.language.lower() == 'chinese' and args.use_softword:
        mhots_config = ConfigDict({'softword': MultiHotConfig(field='softword', use_emb=False)})
    else:
        mhots_config = None

    if args.language.lower() == 'english' and args.char_arch.lower() != 'none':
        char_config = CharConfig(emb_dim=16,
                                 encoder=EncoderConfig(arch=args.char_arch, hid_dim=128, num_layers=1,
                                                       in_drop_rates=(args.drop_rate, 0.0, 0.0)))
        nested_ohots_config = ConfigDict({'char': char_config})
    elif args.language.lower() == 'chinese' and args.use_softlexicon:
        vectors = load_vectors(args.language, 50)
        nested_ohots_config = ConfigDict({'softlexicon': SoftLexiconConfig(vectors=vectors, emb_dim=50, freeze=args.emb_freeze)})
    else:
        nested_ohots_config = None

    if args.use_interm1:
        interm1_config = EncoderConfig(arch=args.enc_arch, hid_dim=args.hid_dim, num_layers=args.num_layers, in_drop_rates=drop_rates)
    else:
        interm1_config = None

    if args.use_interm2:
        interm2_config = EncoderConfig(arch=args.enc_arch, hid_dim=args.hid_dim, num_layers=args.num_layers, in_drop_rates=drop_rates)
    else:
        interm2_config = None

    if args.language.lower() == 'english' and args.use_elmo:
        elmo_config = ELMoConfig(elmo=load_pretrained('elmo'))
    else:
        elmo_config = None

    if args.language.lower() == 'english' and args.use_flair:
        flair_fw_lm, flair_bw_lm = load_pretrained('flair')
        flair_fw_config, flair_bw_config = FlairConfig(flair_lm=flair_fw_lm), FlairConfig(flair_lm=flair_bw_lm)
        if interm2_config is not None:
            interm2_config.in_proj = True
    else:
        flair_fw_config, flair_bw_config = None, None

    if args.bert_arch.lower() != 'none':
        # Cased tokenizer for NER task
        bert_like, tokenizer = load_pretrained(args.bert_arch, args, cased=True)
        bert_like_config = BertLikeConfig(tokenizer=tokenizer, bert_like=bert_like, arch=args.bert_arch,
                                          freeze=False, use_truecase='cased' in os.path.basename(bert_like.name_or_path).split('-'))
    else:
        bert_like_config = None

    return {'ohots': ohots_config,
            'mhots': mhots_config,
            'nested_ohots': nested_ohots_config,
            'intermediate1': interm1_config,
            'elmo': elmo_config,
            'flair_fw': flair_fw_config,
            'flair_bw': flair_bw_config,
            'bert_like': bert_like_config,
            'intermediate2': interm2_config}

def build_ER_config(args: argparse.Namespace):
    drop_rates = (0.0, 0.05, args.drop_rate) if args.use_locked_drop else (args.drop_rate, 0.0, 0.0)

    if args.ck_decoder == 'sequence_tagging':
        decoder_config = SequenceTaggingDecoderConfig(scheme=args.scheme,
                                                      use_crf=args.use_crf,
                                                      fl_gamma=args.fl_gamma,
                                                      sl_epsilon=args.sl_epsilon,
                                                      in_drop_rates=drop_rates)
    elif args.ck_decoder == 'span_classification':
        decoder_config = SpanClassificationDecoderConfig(agg_mode=args.agg_mode,
                                                         fl_gamma=args.fl_gamma,
                                                         sl_epsilon=args.sl_epsilon,
                                                         num_neg_chunks=args.num_neg_chunks,
                                                         max_span_size=args.max_span_size,
                                                         size_emb_dim=args.ck_size_emb_dim,
                                                         in_drop_rates=drop_rates)
    elif args.ck_decoder == 'boundary_selection':
        decoder_config = BoundarySelectionDecoderConfig(use_biaffine=args.use_biaffine,
                                                        affine=EncoderConfig(arch=args.affine_arch, hid_dim=150, num_layers=1, in_drop_rates=(0.4, 0.0, 0.0), hid_drop_rate=0.2),
                                                        fl_gamma=args.fl_gamma,
                                                        sl_epsilon=args.sl_epsilon,
                                                        neg_sampling_rate=args.neg_sampling_rate,
                                                        hard_neg_sampling_rate=args.hard_neg_sampling_rate,
                                                        hard_neg_sampling_size=args.hard_neg_sampling_size,
                                                        sb_epsilon=args.sb_epsilon,
                                                        sb_size=args.sb_size,
                                                        sb_adj_factor=args.sb_adj_factor,
                                                        hid_drop_rates=drop_rates)
    return ExtractorConfig(**collect_IE_assembly_config(args), decoder=decoder_config)

def process_IE_data(train_data, dev_data, test_data, args, config):
    if (config.bert_like is not None and
            ((args.dataset in ('SIGHAN2006', 'yidu_s4k')) or
             (args.dataset in ('conll2003', 'conll2012') and args.doc_level))):
        train_data = segment_uniformly_for_bert_like(train_data, config.bert_like.tokenizer, verbose=args.log_terminal)
        dev_data   = segment_uniformly_for_bert_like(dev_data,   config.bert_like.tokenizer, verbose=args.log_terminal)
        test_data  = segment_uniformly_for_bert_like(test_data,  config.bert_like.tokenizer, verbose=args.log_terminal)

    if args.use_softword or args.use_softlexicon:
        if config.nested_ohots is not None and 'softlexicon' in config.nested_ohots.keys():
            vectors = config.nested_ohots['softlexicon'].vectors
        else:
            vectors = load_vectors(args.language, 50)
        tokenizer = LexiconTokenizer(vectors.itos)
        for data in [train_data, dev_data, test_data]:
            for data_entry in data:
                data_entry['tokens'].build_softwords(tokenizer.tokenize)
                data_entry['tokens'].build_softlexicons(tokenizer.tokenize)

    return train_data, dev_data, test_data

if __name__ == '__main__':
    parser = argparse.ArgumentParser(formatter_class=argparse.ArgumentDefaultsHelpFormatter,
                                     fromfile_prefix_chars='@')
    args = parse_arguments(parser)

    # Use micro-seconds to ensure different timestamps while adopting multiprocessing
    timestamp = datetime.datetime.now().strftime("%Y%m%d-%H%M%S-%f")
    save_path =  f"cache/{args.dataset}-ER/{timestamp}"
    if not os.path.exists(save_path):
        os.makedirs(save_path)

    handlers = [logging.FileHandler(f"{save_path}/training.log")]
    if args.log_terminal:
        handlers.append(logging.StreamHandler(sys.stdout))
    logging.basicConfig(level=logging.INFO,
                        format="[%(asctime)s %(levelname)s] %(message)s",
                        datefmt="%Y-%m-%d %H:%M:%S",
                        handlers=handlers)

    logger = logging.getLogger(__name__)
    logger.info(header_format("Starting", sep='='))
    logger.info(" ".join(sys.argv))
    logger.info(pprint.pformat(args.__dict__))

    logger.info(header_format("Preparing", sep='-'))
    device = auto_device()
    if device.type.startswith('cuda'):
        torch.cuda.set_device(device)
        temp = torch.randn(100).to(device)

    train_data, dev_data, test_data = load_data(args)
    args.language = dataset2language[args.dataset]
    config = build_ER_config(args)
    train_data, dev_data, test_data = process_IE_data(train_data, dev_data, test_data, args, config)

    if not args.train_with_dev:
        train_set = Dataset(train_data, config, training=True)
        train_set.build_vocabs_and_dims(dev_data, test_data)
        dev_set   = Dataset(dev_data,  train_set.config, training=False)
        test_set  = Dataset(test_data, train_set.config, training=False)

        train_loader = torch.utils.data.DataLoader(train_set, batch_size=args.batch_size, shuffle=True,  collate_fn=train_set.collate)
        dev_loader   = torch.utils.data.DataLoader(dev_set,   batch_size=args.batch_size, shuffle=False, collate_fn=dev_set.collate)
    else:
        train_set = Dataset(train_data + dev_data, config, training=True)
        train_set.build_vocabs_and_dims(test_data)
        dev_set   = Dataset([],        train_set.config, training=False)
        test_set  = Dataset(test_data, train_set.config, training=False)

        train_loader = torch.utils.data.DataLoader(train_set, batch_size=args.batch_size, shuffle=True,  collate_fn=train_set.collate)
        dev_loader   = None

    logger.info(train_set.summary)

    logger.info(header_format("Building", sep='-'))
    model = config.instantiate().to(device)
    count_params(model)

    logger.info(header_format("Training", sep='-'))
    trainer = build_trainer(model, device, len(train_loader), args)
    if args.pdb:
        pdb.set_trace()

    if args.profile:
        prof = profile(trainer, train_loader)
        pdb.set_trace()

    torch.save(config, f"{save_path}/{config.name}-config.pth")
    def save_callback(model):
        torch.save(model, f"{save_path}/{config.name}.pth")
    trainer.train_steps(train_loader=train_loader, dev_loader=dev_loader, num_epochs=args.num_epochs,
                        save_callback=save_callback, save_by_loss=False)

    # Save the final version
    # torch.save(model, f"{save_path}/{config.name}.fv.pth")

    logger.info(header_format("Evaluating", sep='-'))
    model = torch.load(f"{save_path}/{config.name}.pth", map_location=device)
    trainer = Trainer(model, device=device)

    logger.info("Evaluating on dev-set")
    evaluate_entity_recognition(trainer, dev_set, batch_size=args.batch_size)
    logger.info("Evaluating on test-set")
    evaluate_entity_recognition(trainer, test_set, batch_size=args.batch_size)

    if args.pipeline:
        # Replace gold chunks with predicted chunks for pipeline
        if args.train_with_dev:
            # Retrieve the original splits
            train_set = Dataset(train_data, train_set.config, training=True)
            dev_set   = Dataset(dev_data,   train_set.config, training=False)

        train_set_chunks_pred = trainer.predict(train_set, batch_size=args.batch_size)
        for ex, chunks_pred in zip(train_data, train_set_chunks_pred):
            ex['chunks'] = ex['chunks'] + [ck for ck in chunks_pred if ck not in ex['chunks']]

        dev_set_chunks_pred = trainer.predict(dev_set, batch_size=args.batch_size)
        for ex, chunks_pred in zip(dev_data, dev_set_chunks_pred):
            ex['chunks'] = chunks_pred

        test_set_chunks_pred = trainer.predict(test_set, batch_size=args.batch_size)
        for ex, chunks_pred in zip(test_data, test_set_chunks_pred):
            ex['chunks'] = chunks_pred

        torch.save((train_data, dev_data, test_data), f"{save_path}/data-with-chunks-pred.pth")

    logger.info(" ".join(sys.argv))
    logger.info(pprint.pformat(args.__dict__))
    logger.info(header_format("Ending", sep='='))

max_stars_repo_name: syuoni/eznlp

max_issues_repo_name: syuoni/eznlp

max_forks_repo_name: syuoni/eznlp