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

from cupy import (
    arange,
    array,
    asarray,
    atleast_1d,
    intc,
    integer,
    isscalar,
    issubdtype,
    take,
    unique,
    where,
)
from bisect import bisect_left

def next_fast_len(target):
    """
    Find the next fast size of input data to `fft`, for zero-padding, etc.

    SciPy's FFTPACK has efficient functions for radix {2, 3, 4, 5}, so this
    returns the next composite of the prime factors 2, 3, and 5 which is
    greater than or equal to `target`. (These are also known as 5-smooth
    numbers, regular numbers, or Hamming numbers.)

    Parameters
    ----------
    target : int
        Length to start searching from.  Must be a positive integer.

    Returns
    -------
    out : int
        The first 5-smooth number greater than or equal to `target`.

    Notes
    -----
    .. versionadded:: 0.18.0

    Examples
    --------
    On a particular machine, an FFT of prime length takes 133 ms:

    >>> from scipy import fftpack
    >>> min_len = 10007  # prime length is worst case for speed
    >>> a = np.random.randn(min_len)
    >>> b = fftpack.fft(a)

    Zero-padding to the next 5-smooth length reduces computation time to
    211 us, a speedup of 630 times:

    >>> fftpack.helper.next_fast_len(min_len)
    10125
    >>> b = fftpack.fft(a, 10125)

    Rounding up to the next power of 2 is not optimal, taking 367 us to
    compute, 1.7 times as long as the 5-smooth size:

    >>> b = fftpack.fft(a, 16384)

    """
    hams = (
        8,
        9,
        10,
        12,
        15,
        16,
        18,
        20,
        24,
        25,
        27,
        30,
        32,
        36,
        40,
        45,
        48,
        50,
        54,
        60,
        64,
        72,
        75,
        80,
        81,
        90,
        96,
        100,
        108,
        120,
        125,
        128,
        135,
        144,
        150,
        160,
        162,
        180,
        192,
        200,
        216,
        225,
        240,
        243,
        250,
        256,
        270,
        288,
        300,
        320,
        324,
        360,
        375,
        384,
        400,
        405,
        432,
        450,
        480,
        486,
        500,
        512,
        540,
        576,
        600,
        625,
        640,
        648,
        675,
        720,
        729,
        750,
        768,
        800,
        810,
        864,
        900,
        960,
        972,
        1000,
        1024,
        1080,
        1125,
        1152,
        1200,
        1215,
        1250,
        1280,
        1296,
        1350,
        1440,
        1458,
        1500,
        1536,
        1600,
        1620,
        1728,
        1800,
        1875,
        1920,
        1944,
        2000,
        2025,
        2048,
        2160,
        2187,
        2250,
        2304,
        2400,
        2430,
        2500,
        2560,
        2592,
        2700,
        2880,
        2916,
        3000,
        3072,
        3125,
        3200,
        3240,
        3375,
        3456,
        3600,
        3645,
        3750,
        3840,
        3888,
        4000,
        4050,
        4096,
        4320,
        4374,
        4500,
        4608,
        4800,
        4860,
        5000,
        5120,
        5184,
        5400,
        5625,
        5760,
        5832,
        6000,
        6075,
        6144,
        6250,
        6400,
        6480,
        6561,
        6750,
        6912,
        7200,
        7290,
        7500,
        7680,
        7776,
        8000,
        8100,
        8192,
        8640,
        8748,
        9000,
        9216,
        9375,
        9600,
        9720,
        10000,
    )

    target = int(target)

    if target <= 6:
        return target

    # Quickly check if it's already a power of 2
    if not (target & (target - 1)):
        return target

    # Get result quickly for small sizes, since FFT itself is similarly fast.
    if target <= hams[-1]:
        return hams[bisect_left(hams, target)]

    match = float("inf")  # Anything found will be smaller
    p5 = 1
    while p5 < target:
        p35 = p5
        while p35 < target:
            # Ceiling integer division, avoiding conversion to float
            # (quotient = ceil(target / p35))
            quotient = -(-target // p35)

            # Quickly find next power of 2 >= quotient
            p2 = 2 ** ((quotient - 1).bit_length())

            N = p2 * p35
            if N == target:
                return N
            elif N < match:
                match = N
            p35 *= 3
            if p35 == target:
                return p35
        if p35 < match:
            match = p35
        p5 *= 5
        if p5 == target:
            return p5
    if p5 < match:
        match = p5
    return match

def _init_nd_shape_and_axes(x, shape, axes):
    """Handle shape and axes arguments for n-dimensional transforms.

    Returns the shape and axes in a standard form, taking into account negative
    values and checking for various potential errors.

    Parameters
    ----------
    x : array_like
        The input array.
    shape : int or array_like of ints or None
        The shape of the result.  If both `shape` and `axes` (see below) are
        None, `shape` is ``x.shape``; if `shape` is None but `axes` is
        not None, then `shape` is ``scipy.take(x.shape, axes, axis=0)``.
        If `shape` is -1, the size of the corresponding dimension of `x` is
        used.
    axes : int or array_like of ints or None
        Axes along which the calculation is computed.
        The default is over all axes.
        Negative indices are automatically converted to their positive
        counterpart.

    Returns
    -------
    shape : array
        The shape of the result. It is a 1D integer array.
    axes : array
        The shape of the result. It is a 1D integer array.

    """
    x = asarray(x)
    noshape = shape is None
    noaxes = axes is None

    if noaxes:
        axes = arange(x.ndim, dtype=intc)
    else:
        axes = atleast_1d(axes)

    if axes.size == 0:
        axes = axes.astype(intc)

    if not axes.ndim == 1:
        raise ValueError("when given, axes values must be a scalar or vector")
    if not issubdtype(axes.dtype, integer):
        raise ValueError("when given, axes values must be integers")

    axes = where(axes < 0, axes + x.ndim, axes)

    if axes.size != 0 and (axes.max() >= x.ndim or axes.min() < 0):
        raise ValueError("axes exceeds dimensionality of input")
    if axes.size != 0 and unique(axes).shape != axes.shape:
        raise ValueError("all axes must be unique")

    if not noshape:
        shape = atleast_1d(shape)
    elif isscalar(x):
        shape = array([], dtype=intc)
    elif noaxes:
        shape = array(x.shape, dtype=intc)
    else:
        shape = take(x.shape, axes)

    if shape.size == 0:
        shape = shape.astype(intc)

    if shape.ndim != 1:
        raise ValueError("when given, shape values must be a scalar or vector")
    if not issubdtype(shape.dtype, integer):
        raise ValueError("when given, shape values must be integers")
    if axes.shape != shape.shape:
        raise ValueError(
            "when given, axes and shape arguments"
            " have to be of the same length"
        )

    shape = where(shape == -1, array(x.shape)[axes], shape)

    if shape.size != 0 and (shape < 1).any():
        raise ValueError(
            "invalid number of data points ({0}) specified".format(shape)
        )

    return shape, axes

def _init_nd_shape_and_axes_sorted(x, shape, axes):
    """Handle and sort shape and axes arguments for n-dimensional transforms.

    This is identical to `_init_nd_shape_and_axes`, except the axes are
    returned in sorted order and the shape is reordered to match.

    Parameters
    ----------
    x : array_like
        The input array.
    shape : int or array_like of ints or None
        The shape of the result.  If both `shape` and `axes` (see below) are
        None, `shape` is ``x.shape``; if `shape` is None but `axes` is
        not None, then `shape` is ``scipy.take(x.shape, axes, axis=0)``.
        If `shape` is -1, the size of the corresponding dimension of `x` is
        used.
    axes : int or array_like of ints or None
        Axes along which the calculation is computed.
        The default is over all axes.
        Negative indices are automatically converted to their positive
        counterpart.

    Returns
    -------
    shape : array
        The shape of the result. It is a 1D integer array.
    axes : array
        The shape of the result. It is a 1D integer array.

    """
    x = asarray(x)
    noaxes = axes is None
    shape, axes = _init_nd_shape_and_axes(x, shape, axes)

    if not noaxes:
        shape = shape[axes.argsort()]
        axes.sort()

    return shape, axes

reason: from scipy,from cupy


text:
import numpy as np

def solidangle(twotheta, dtwotheta, sd, dsd, pixelsize, dpixelsize):
    """Solid-angle correction for two-dimensional SAS images with error propagation

    Inputs:
        twotheta: matrix of two-theta values
        dtwotheta: matrix of absolute error of two-theta values
        sampletodetectordistance: sample-to-detector distance
        dsampletodetectordistance: absolute error of sample-to-detector distance

    Outputs two matrices of the same shape as twotheta. The scattering intensity
        matrix should be multiplied by the first one. The second one is the propagated
        error of the first one.
    """
    sac = sd ** 2 / np.cos(twotheta) ** 3 / pixelsize ** 2
    return (sac,
            sac * (4 * dsd ** 2 / sd ** 2 + 4 * dpixelsize ** 2 / pixelsize ** 2 + 9 * dtwotheta ** 2 * np.tan(
                twotheta) ** 2) ** 0.5)

def _angledependentabsorption_value(twotheta: np.ndarray, transmission: float):
    cor = np.ones_like(twotheta)
    if transmission == 1:
        return cor
    mud = -np.log(transmission)

    cor[twotheta > 0] = transmission * mud * (1 - 1 / np.cos(twotheta[twotheta > 0])) / (
            np.exp(-mud / np.cos(twotheta[twotheta > 0])) - transmission)
    return cor

def _angledependentabsorption_error(twotheta, dtwotheta, transmission, dtransmission):
    # calculated using sympy
    costth = np.cos(twotheta)
    sintth = np.sin(twotheta)
    lntrans = np.log(transmission)
    exp1 = np.exp(lntrans / costth)
    return ((transmission * costth - exp1 *
             lntrans * costth + exp1
             * lntrans - exp1 * costth) ** 2
            * (transmission ** 2 * dtwotheta ** 2 * lntrans ** 2 * sintth ** 2
               + dtransmission ** 2 * sintth ** 4 - 3 * dtransmission ** 2 * sintth ** 2
               - 2 * dtransmission ** 2 * costth ** 3 + 2 * dtransmission ** 2) /
            (transmission - exp1) ** 4) ** 0.5 * np.abs(costth) ** (-3.0)

def angledependentabsorption(twotheta, dtwotheta, transmission, dtransmission):
    """Correction for angle-dependent absorption of the sample with error propagation

    Inputs:
        twotheta: matrix of two-theta values
        dtwotheta: matrix of absolute error of two-theta values
        transmission: the transmission of the sample (I_after/I_before, or
            exp(-mu*d))
        dtransmission: the absolute error of the transmission of the sample

    Two matrices are returned: the first one is the correction (intensity matrix
        should be multiplied by it), the second is its absolute error.
    """
    # error propagation formula calculated using sympy
    return (_angledependentabsorption_value(twotheta, transmission),
            _angledependentabsorption_error(twotheta, dtwotheta, transmission, dtransmission))

def angledependentairtransmission(twotheta, dtwotheta,
                                  pressure, sampletodetectordistance,
                                  dsampletodetectordistance, mu0_air=1 / 883.49, dmu0_air=0):
    """Correction for the angle dependent absorption of air in the scattered
    beam path, with error propagation

    Inputs:
            twotheta: matrix of two-theta values
            dtwotheta: absolute error matrix of two-theta
            pressure: the air pressure in mbar
            sampletodetectordistance: sample-to-detector distance
            dsampletodetectordistance: error of the sample-to-detector distance
            mu0_air: the linear absorption coefficient of air at 1000 mbars, in 1/mm units.
                At 8 keV this is 1/883.49 1/mm
            dmu0_air: error of the linear absorption coefficient of air at 1000 mbars in 1/mm units

    1/mu_air and sampletodetectordistance should have the same dimension

    The scattering intensity matrix should be multiplied by the resulting
    correction matrix."""
    if (pressure <= 0) or (not np.isfinite(pressure)):
        return np.ones_like(twotheta), np.zeros_like(dtwotheta)
    costth = np.cos(twotheta)
    sintth = np.sin(twotheta)
    mu_air = mu0_air / 1000 * pressure
    dmu_air = dmu0_air / 1000 * pressure
    return (np.exp(mu_air * sampletodetectordistance / costth),
            np.sqrt(dmu_air ** 2 * sampletodetectordistance ** 2 *
                    np.exp(2 * mu_air * sampletodetectordistance / costth)
                    / costth ** 2 + dsampletodetectordistance ** 2 *
                    mu_air ** 2 * np.exp(2 * mu_air * sampletodetectordistance /
                                         costth) /
                    costth ** 2 + dtwotheta ** 2 * mu_air ** 2 *
                    sampletodetectordistance ** 2 *
                    np.exp(2 * mu_air * sampletodetectordistance / costth)
                    * sintth ** 2 / costth ** 4)
            )

ext: py

max_stars_repo_name: awacha/cct

max_issues_repo_name: awacha/cct

max_forks_repo_name: awacha/cct

reason: import numpy


text:
# ニューラルネットワーク
# 一番簡単なBP
# XOR

import numpy as np

EPSILON = 1.0
ETA = 0.1
TIME = 10000

def sigmoid(x):
    return 1/(1+np.exp(-EPSILON*x))

inputs = np.array([[0,0],
         [0,1],
         [1,0],
         [1,1]])
teach = np.array([0,1,1,0])

inputs = inputs.astype(np.float32)
teach = teach.astype(np.float32)

weight1 = np.array([
    [-0.1280102 , -0.94814754],
    [ 0.09932496, -0.12935521]])
weight2 = np.array([[-0.1592644 , -0.33933036]])
offset1 = np.array([-0.59070273,  0.23854193])
offset2 = np.array([-0.40069065])

for t in range(TIME):
    error = 0
    out = []

    w1 = np.zeros((2,2))
    w2 = np.zeros((1,2))
    o1 = np.zeros(2)
    o2 = np.zeros(1)

    for p in range(len(inputs)):
        # 前向き計算
        out1 = sigmoid(np.dot(weight1,inputs[p])+offset1)
        out2 = sigmoid(np.dot(weight2, out1)+offset2)
        out.append(out2)
        error += (out2-teach[p])**2

        # BP
        delta2 = (out2-teach[p])*EPSILON*out2*(1.0-out2)
        w2 -= ETA*delta2*out1
        o2 -= ETA*delta2

        delta1 = EPSILON*out1*(1.0-out1)*delta2*weight2
        w1 -= ETA*delta1*inputs[p]
        o1 -= ETA*delta1[0]

    weight2 += w2
    weight1 += w1
    offset2 += o2
    offset1 += o1

    print(error)

print('output')
print(out)


text:
import json
from dataclasses import dataclass
from datetime import datetime
from random import shuffle, choice
from typing import Dict, List, Iterable, Mapping, Tuple, Optional

import math
import numpy
import trueskill
import itertools
from pathlib import Path
from rlbot.parsing.bot_config_bundle import BotConfigBundle, get_bot_config_bundle

from bots import BotID, fmt_bot_name
from leaguesettings import LeagueSettings
from match import MatchDetails
from paths import LeagueDir, PackageFiles
from ranking_system import RankingSystem
from trueskill import Rating

# Minimum required TrueSkill match quality. Can't be higher than 0.44
MIN_REQ_FAIRNESS = 0.3

class TicketSystem:
    def __init__(self):
        self.tickets: Dict[BotID, float] = {}
        self.new_bot_ticket_count = 4.0
        self.session_game_counts: Dict[BotID, int] = {}

        # Decrease this number toward 1.0 if you want to prioritize a balanced number of games played.
        # Increase it if you want more randomness, and priority for bots who haven't played recently.
        self.ticket_increase_rate = 1.2

        # Increase this number to more heavily prioritize bots who have not played many games yet
        # during the current session. Can be anything >= 0.
        self.game_catchup_boost = 1.0

    def ensure(self, bots: Iterable[BotID]):
        """
        Ensure that all bots in the given list has tickets in the ticket system.
        """
        for bot in bots:
            if bot not in self.tickets:
                # Give new bots some tickets right away
                self.tickets[bot] = self.new_bot_ticket_count
            if bot not in self.session_game_counts:
                self.session_game_counts[bot] = 0

    def get_ensured(self, bot: BotID) -> float:
        """
        Returns the number of tickets owned by the given bot. The bot is added with the default
        number of tickets, if they are not in the system yet.
        """
        if bot not in self.tickets:
            self.tickets[bot] = self.new_bot_ticket_count
        return self.tickets[bot]

    def get(self, bot: BotID) -> Optional[float]:
        """
        Returns the number of tickets owned by the given bot or None of the bot is not in the system.
        """
        return self.tickets.get(bot)

    def set(self, bot: BotID, tickets: float):
        """
        Set the number of tickets for the given bot. The update is not persistent until `save` is called.
        """
        self.tickets[bot] = tickets

    def total(self) -> float:
        """
        Returns the total number of tickets in the ticket system.
        """
        return sum(self.tickets.values())

    def pick_bots(self, bots: Iterable[BotID]) -> List[BotID]:
        """
        Picks 6 unique bots based on their number of tickets in the ticket system
        """
        self.ensure(bots)

        # We don't use self.total() since it can be the case, that not all bots appear in `bots`
        bot_tickets = [self.get_ensured(bot_id) for bot_id in bots]
        total = sum(bot_tickets)
        prop = [tickets / total for tickets in bot_tickets]
        picked = list(numpy.random.choice(list(bots), 6, p=prop, replace=False))

        return picked

    def choose(self, chosen_bots: Iterable[BotID], all_bots: Iterable[BotID]):
        """
        Choose the list of given bots, which will reset their number of tickets and double every else's.
        """
        max_game_count = max(self.session_game_counts.values())
        for bot in all_bots:
            if bot in chosen_bots:
                # Reset their tickets
                self.tickets[bot] = 1.0
                self.session_game_counts[bot] += 1
            else:
                # Increase their tickets

                # Tickets increase faster if the bot is lagging behind on the number of games played.
                games_deficit = max_game_count - self.session_game_counts[bot]

                # Tickets also multiply a little even if the bot has played more games than any other.
                self.tickets[bot] *= (self.ticket_increase_rate + games_deficit * self.game_catchup_boost)

    def save(self, ld: LeagueDir, time_stamp: str):
        with open(ld.tickets / f"{time_stamp}_tickets.json", 'w') as f:
            json.dump(self.tickets, f, sort_keys=True)

    @staticmethod
    def load(ld: LeagueDir) -> 'TicketSystem':
        ticket_sys = TicketSystem()
        if any(ld.tickets.iterdir()):
            # Assume last tickets file is the newest, since they are prefixed with a time stamp
            with open(list(ld.tickets.iterdir())[-1]) as f:
                ticket_sys.tickets = json.load(f)

        settings = LeagueSettings.load(ld)
        ticket_sys.new_bot_ticket_count = settings.new_bot_ticket_count
        ticket_sys.ticket_increase_rate = settings.ticket_increase_rate
        ticket_sys.game_catchup_boost = settings.game_catchup_boost

        matches_in_session = MatchDetails.latest(ld, settings.last_summary)
        for match in matches_in_session:
            bots = match.blue + match.orange
            ticket_sys.ensure(bots)
            for bot_id in bots:
                ticket_sys.session_game_counts[bot_id] += 1

        return ticket_sys

    @staticmethod
    def read(path: Path, settings: LeagueSettings) -> 'TicketSystem':
        ticket_sys = TicketSystem()
        with open(path) as f:
            ticket_sys.tickets = json.load(f)
            ticket_sys.new_bot_ticket_count = settings.new_bot_ticket_count
            ticket_sys.ticket_increase_rate = settings.ticket_increase_rate
            ticket_sys.game_catchup_boost = settings.game_catchup_boost
            return ticket_sys

    @staticmethod
    def all(ld: LeagueDir, settings: LeagueSettings):
        first = TicketSystem()
        first.new_bot_ticket_count = settings.new_bot_ticket_count
        first.ticket_increase_rate = settings.ticket_increase_rate
        first.game_catchup_boost = settings.game_catchup_boost
        return [first] + [TicketSystem.read(path, settings) for path in list(ld.tickets.iterdir())]

    @staticmethod
    def undo(ld: LeagueDir):
        """
        Remove latest tickets file
        """
        if any(ld.tickets.iterdir()):
            # Assume last tickets file is the newest, since they are prefixed with a time stamp
            list(ld.tickets.iterdir())[-1].unlink()  # Remove file
        else:
            print("No tickets to undo.")

@dataclass
class Candidate:
    bot_id: BotID
    rating: Rating

# We need this here until Trueskill fixes their implementation
def pdf(x, mu: float = 0, sigma: float = 1):
    """Probability density function"""
    return (1 / (math.sqrt(2 * math.pi) * abs(sigma)) *
            math.exp(-(((x - mu) / abs(sigma)) ** 2 / 2)))

class MatchMaker:
    @staticmethod
    def make_next(bots: Mapping[BotID, BotConfigBundle], rank_sys: RankingSystem,
                  ticket_sys: TicketSystem) -> MatchDetails:
        """
        Make the next match to play. This will use to TicketSystem and the RankingSystem to find
        a fair match between some bots that haven't played for a while. It is assumed that the match
        is guaranteed to finish (since the TicketSystem is updated).
        """

        time_stamp = make_timestamp()
        blue, orange = MatchMaker.decide_on_players_3(bots.keys(), rank_sys, ticket_sys)
        name = "_".join([time_stamp] + blue + ["vs"] + orange)
        map = choice([
            "ChampionsField",
            "DFHStadium",
            "NeoTokyo",
            "UrbanCentral",
            "BeckwithPark",
        ])
        return MatchDetails(time_stamp, name, blue, orange, map)

    @staticmethod
    def decide_on_players(bot_ids: Iterable[BotID], rank_sys: RankingSystem,
                          ticket_sys: TicketSystem) -> Tuple[List[BotID], List[BotID]]:
        """
        Find two balanced teams. The TicketSystem and the RankingSystem to find
        a fair match up between some bots that haven't played for a while.
        """
        limit = 400

        tries_left = limit
        while tries_left > 0:
            tries_left -= 1

            # Pick some bots that haven't played for a while
            picked = ticket_sys.pick_bots(bot_ids)
            shuffle(picked)
            ratings = [rank_sys.get(bot) for bot in picked]

            blue = tuple(ratings[0:3])
            orange = tuple(ratings[3:6])

            # Is this a fair match?
            required_fairness = min(tries_left / limit, MIN_REQ_FAIRNESS)
            if trueskill.quality([blue, orange]) >= required_fairness:
                tickets_consumed = sum([ticket_sys.get_ensured(b) for b in picked])
                print(f"Match: {picked[0:3]} vs {picked[3:6]}\nMatch quality: {trueskill.quality([blue, orange])}  Tickets consumed: {tickets_consumed}")
                ticket_sys.choose(picked, bot_ids)
                return picked[0:3], picked[3:6]

        raise Exception("Failed to find a fair match")

    @staticmethod
    def decide_on_players_2(bot_ids: Iterable[BotID], rank_sys: RankingSystem,
                            ticket_sys: TicketSystem) -> Tuple[List[BotID], List[BotID]]:
        """
        Find two balanced teams. The TicketSystem and the RankingSystem to find
        a fair match up between some bots that haven't played for a while.
        """

        # Composing a team of the best player + the worst two players will likely yield a balanced match (0, 4, 5).
        # These represent a few arrangements like that which seem reasonable to try, they will be checked against
        # the trueskill system.
        likely_balances = [(0, 4, 5), (0, 3, 5), (0, 2, 5), (0, 3, 4)]

        # Experimental average quality based on limit:
        # 1000: 0.4615
        # 400:  0.460
        # 100:  0.457
        # 10:   0.448
        num_bot_groups_to_test = 400

        # How much we value the tightness of rating distribution in a given match.
        # A higher number will yield matches with similarly skilled bots, but potentially lower probability of a draw.
        tightness_weight = 1.0

        tries_left = num_bot_groups_to_test
        best_quality_found = 0
        best_score_found = 0
        best_match = None
        chosen_balance = None

        while tries_left > 0:
            tries_left -= 1

            # Pick some bots that haven't played for a while
            picked = ticket_sys.pick_bots(bot_ids)
            candidates = [Candidate(bot, rank_sys.get(bot)) for bot in picked]
            candidates.sort(key=lambda c: float(c.rating), reverse=True)
            tightness = 1 / (numpy.std([float(c.rating) for c in candidates]) + 1)

            for balance in likely_balances:
                blue_candidates = candidates[balance[0]], candidates[balance[1]], candidates[balance[2]]
                orange_candidates = [c for c in candidates if c not in blue_candidates]
                quality = trueskill.quality([[c.rating for c in blue_candidates], [c.rating for c in orange_candidates]])
                score = quality + tightness * tightness_weight
                if score > best_score_found:
                    best_score_found = score
                    best_quality_found = quality
                    best_match = (blue_candidates, orange_candidates)
                    chosen_balance = balance

        blue_ids = [c.bot_id for c in best_match[0]]
        orange_ids = [c.bot_id for c in best_match[1]]
        tickets_consumed = sum([ticket_sys.get_ensured(b) for b in blue_ids + orange_ids])
        print(f"Match: {blue_ids} vs {orange_ids}\nMatch quality: {best_quality_found}  score: {best_score_found}  "
              f"Rank pattern: {chosen_balance}")
        ticket_sys.choose(blue_ids + orange_ids, bot_ids)
        return blue_ids, orange_ids

    @staticmethod
    def decide_on_players_3(bot_ids: Iterable[BotID], rank_sys: RankingSystem,
                            ticket_sys: TicketSystem) -> Tuple[List[BotID], List[BotID]]:
        """
        Find two balanced teams. The TicketSystem and the RankingSystem to find
        a fair match up between some bots that haven't played for a while.
        """
        # Higher ticket strength produces a more uniform distribution of matches played, adjust by increments of 0.1
        TICKET_STRENGTH = 1
        # Higher MMR tolerance allows accurately rated bots to play in more "distant" MMR matches, adjust by increments of 1
        MMR_TOLERANCE = 4
        # Max attempts to build match of quality >= MIN_QUALITY
        MAX_ITERATIONS = 20
        MIN_QUALITY = 0.4

        rank_sys.ensure_all(bot_ids)
        ticket_sys.ensure(bot_ids)

        best_quality = 0
        best_match = None

        for i in range(MAX_ITERATIONS):
            # Get Leader Bot (choose randomly between bots with highest tickets)
            max_tickets = max(ticket_sys.tickets.values())
            possible_leaders = [bot_id for bot_id, tickets in ticket_sys.tickets.items() if tickets == max_tickets]
            leader = numpy.random.choice(possible_leaders)

            # Get MU for Leader bot, that will be the match mmr
            match_mmr = rank_sys.get(leader).mu

            # Score all bots based on probability to perform at target mmr, scaled by amount of tickets
            candidates = [Candidate(bot_id, rank_sys.get(bot_id)) for bot_id in bot_ids if bot_id != leader]
            scores = []

            for c in candidates:
                # Calculate probability to perform at desired mmr
                performance_prob = pdf(match_mmr, mu=c.rating.mu, sigma=math.sqrt(c.rating.sigma**2 + MMR_TOLERANCE**2))

                # Calculate weighting factor based on tickets
                tickets = ticket_sys.get(c.bot_id)
                tickets_weight = tickets ** TICKET_STRENGTH

                # Calculate candidate score
                scores.append(performance_prob * tickets_weight)

            # Pick 5 bots randomly based on their score
            probs = numpy.asarray(scores) / sum(scores)
            players = list(numpy.random.choice(candidates, size=5, p=probs, replace=False))
            players.append(Candidate(leader, rank_sys.get(leader)))

            # Get the highest quality match with the 6 chosen bots
            combinations = list(itertools.combinations(players, 3))
            possible_matches = len(combinations) // 2
            blue_combs = combinations[:possible_matches]
            orange_combs = combinations[:possible_matches-1:-1]

            for i in range(possible_matches):
                blue_team = blue_combs[i]
                orange_team = orange_combs[i]
                quality = trueskill.quality([[c.rating for c in blue_team], [c.rating for c in orange_team]])
                if quality > best_quality:
                    best_quality = quality
                    best_match = (blue_team, orange_team)

            if best_quality >= MIN_QUALITY:
                break

        # We sort by get_mmr() because it considers sigma
        blue_ids = sorted([c.bot_id for c in best_match[0]], key=lambda id: rank_sys.get_mmr(id), reverse=True)
        orange_ids = sorted([c.bot_id for c in best_match[1]], key=lambda id: rank_sys.get_mmr(id), reverse=True)

        tickets_consumed = sum([ticket_sys.get_ensured(b) for b in blue_ids + orange_ids])
        print(f"Match: {blue_ids} vs {orange_ids}\nMatch quality: {best_quality}  Tickets consumed: {tickets_consumed}")
        ticket_sys.choose(blue_ids + orange_ids, bot_ids)
        return blue_ids, orange_ids

    @staticmethod
    def make_test_match(bot_id: BotID) -> MatchDetails:
        allstar_config = get_bot_config_bundle(PackageFiles.psyonix_allstar)
        allstar_id = fmt_bot_name(allstar_config.name)
        team = [bot_id, allstar_id, allstar_id]
        return MatchDetails("", f"test_{bot_id}", team, team, "ChampionsField")

def make_timestamp() -> str:
    return datetime.now().strftime("%Y%m%d%H%M%S")


text:
# coding: utf-8

# In[1]:

import os
import sys
import datetime as dt

sys.path.append(r"..\..")
import pyphenocam

import matplotlib as mpl
import matplotlib.pyplot as plt

sys.path.append(r"J:\Projects\NCCSC\phenocam\Tools\DaymetPy\daymetpy")
import daymetpy

import rasterio

# In[2]:

base_dname = r"J:\Projects\NCCSC\phenocam\DerivedData\nationalelkrefuge"
site_name = "nationalelkrefuge"

# In[3]:

phenosite = pyphenocam.dataaccess.get_site(site_name)

# In[4]:

landsat_dname = os.path.join(base_dname, 'Landsat', 'SceneSubset')
ndvi_fnames = [f for f in os.listdir(landsat_dname ) if 'ndvi' in f]
gcc_fnames = [f for f in os.listdir(landsat_dname ) if 'gcc' in f]
cloud_fnames = [f for f in os.listdir(landsat_dname ) if 'cloud' in f]

scene_data = {}
for scene in ndvi_fnames+gcc_fnames:
    year = scene[9:13]
    jday = scene[13:16]
    date = dt.datetime(int(year), 1, 1, 12) + dt.timedelta(int(jday))
    which = scene.split('_')[1].split('.')[0]
    if date not in scene_data:
        scene_data[date] = {}
    data = rasterio.open(os.path.join(landsat_dname, scene)).read(1, masked=True)
    cloud_data = rasterio.open(os.path.join(landsat_dname, scene).replace(which, 'cloud')).read(1)
    data.mask += cloud_data==4
    scene_data[date][which] = data

scene_dates= list(scene_data.keys())
scene_dates.sort()

# In[6]:

get_ipython().magic('matplotlib inline')
import matplotlib as mpl
import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1 import ImageGrid
import numpy as np

fig = plt.figure(1, (12., 12.))
grid = ImageGrid(fig, 111,  # similar to subplot(111)
                 nrows_ncols=(6, 6),  # creates 2x2 grid of axes
                 axes_pad=0.1,  # pad between axes in inch.
                 share_all=True,
                    cbar_location="bottom",
                    cbar_mode="single",)
for ax in grid:
    pyphenocam.plotting.format_photo_axes(ax)

for i, scene_date in enumerate(scene_dates):
    im =grid[i].imshow(scene_data[scene_date]['ndvi'], vmin=0., vmax=0.7, cmap=mpl.cm.summer_r)
#     pyphenocam.plotting.add_inner_title(grid[i], scene_date.strftime("%Y-%m-%d"), 9, size=8)
    grid[i].set_title(scene_date.strftime("%Y-%m-%d"), fontsize=10)

grid.cbar_axes[0].colorbar(im)

fig.suptitle("Landsat8 Luo-Trishchenko-Khlopenkov clouds masked", fontsize=20, y=.95)
# plt.tight_layout()

fig.savefig(os.path.join(r"J:\Projects\NCCSC\phenocam\Doc\Presentation", "landsatscenes.jpg"), dpi=270)

# In[ ]:

get_ipython().magic('matplotlib notebook')
import matplotlib as mpl
import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1 import ImageGrid
import numpy as np

fig = plt.figure(1, (12., 12.))
grid = ImageGrid(fig, 111,  # similar to subplot(111)
                 nrows_ncols=(5, 5),  # creates 2x2 grid of axes
                 axes_pad=0.2,  # pad between axes in inch.
                 share_all=True,
                    cbar_location="bottom",
                    cbar_mode="single",)

for i, scene_date in enumerate(scene_dates):
    im =grid[i].imshow(scene_data[scene_date], vmin=0., vmax=.5, cmap=mpl.cm.summer_r)
    grid[i].set_title(scene_date.strftime("%Y-%m-%d"), fontsize=10)
    pyphenocam.plotting.format_photo_axes(grid[i])
grid.cbar_axes[0].colorbar(im)

# In[11]:

base_dname = r"J:\Projects\NCCSC\phenocam\DerivedData\nationalelkrefuge"
site_name = "nationalelkrefuge"
site = pyphenocam.dataaccess.get_site(site_name)

# In[15]:

get_ipython().magic('matplotlib inline')
landsat_fishnet_fname = os.path.join(base_dname, "ArcScene", "landsat_fishnet.bmp")
landsat_index_fname = os.path.join(base_dname, "ArcScene", "landsat_subset_index.bmp")

phenosite = pyphenocam.dataaccess.get_site(site_name)

closest_date = dt.datetime(2016, 5, 29, 12, 30)
print(closest_date)
sample_photo_fname = phenosite.get_closest_fname(closest_date)

local_fname = phenosite.get_local_image_fname(sample_photo_fname)
local_fname_ir = phenosite.get_local_image_fname(sample_photo_fname, IR=True)
sample_image = phenosite.get_local_image(sample_photo_fname)
sample_image_ir = phenosite.get_local_image(sample_photo_fname, IR=True)

fig, ax = plt.subplots(1, figsize=(12,5))
ax.imshow(sample_image)
pyphenocam.plotting.format_photo_axes(ax)

# In[6]:

fig = plt.figure(1, (12., 12.))
grid = ImageGrid(fig, 111,  # similar to subplot(111)
                 nrows_ncols=(4, 4),  # creates 2x2 grid of axes
                 axes_pad=0.2,  # pad between axes in inch.
                 share_all=True,
                    cbar_location="bottom",
                    cbar_mode="single",)

for i, scene_date in enumerate(scene_dates):
    print(scene_date, end=' ')
    closest_photo_fname = phenosite.get_closest_fname(scene_date)
    print(closest_photo_fname)
    im =grid[i].imshow(phenosite.get_local_image(closest_photo_fname))
    pyphenocam.plotting.format_photo_axes(grid[i])
grid.cbar_axes[0].colorbar(im)

# In[ ]:

max_stars_repo_name: talbertc-usgs/pyphenocam

max_issues_repo_name: talbertc-usgs/pyphenocam

max_forks_repo_name: talbertc-usgs/pyphenocam


text:
# This file is part of the pyMOR project (http://www.pymor.org).
# Copyright 2013-2016 pyMOR developers and contributors. All rights reserved.
# License: BSD 2-Clause License (http://opensource.org/licenses/BSD-2-Clause)

import numpy as np

from pymor.algorithms.image import estimate_image_hierarchical
from pymor.core.interfaces import ImmutableInterface
from pymor.core.exceptions import ImageCollectionError
from pymor.core.logger import getLogger
from pymor.operators.basic import OperatorBase
from pymor.operators.constructions import induced_norm
from pymor.reductors.basic import GenericRBReconstructor
from pymor.vectorarrays.numpy import NumpyVectorSpace

def reduce_residual(operator, rhs=None, RB=None, rhs_is_functional=True, product=None, extends=None):
    """Generic reduced basis residual reductor.

    Given an operator and a right-hand side, the residual is given by::

        residual.apply(U, mu) == operator.apply(U, mu) - rhs.as_vector(mu)

    When the rhs is a functional we are usually interested in the Riesz representative
    of the residual::

        residual.apply(U, mu)
            == product.apply_inverse(operator.apply(U, mu) - rhs.as_vector(mu))

    Given a basis `RB` of a subspace of the source space of `operator`, this method
    uses :func:`~pymor.algorithms.image.estimate_image_hierarchical` to determine
    a low-dimensional subspace containing the image of the subspace under
    `residual` (resp. `riesz_residual`), computes an orthonormal basis
    `residual_range` for this range space and then returns the Petrov-Galerkin projection ::

        projected_residual
            === residual.projected(range_basis=residual_range, source_basis=RB)

    of the residual operator. Given an reduced basis coefficient vector `u`, w.r.t.
    `RB`, the (dual) norm of the residual can then be computed as ::

        projected_residual.apply(u, mu).l2_norm()

    Moreover, a `reconstructor` is provided such that ::

        reconstructor.reconstruct(projected_residual.apply(u, mu))
            == residual.apply(RB.lincomb(u), mu)

    Parameters
    ----------
    operator
        See definition of `residual`.
    rhs
        See definition of `residual`. If `None`, zero right-hand side is assumed.
    rhs_is_functional
        Set this to `True` when `rhs` is a |Functional|.
    RB
        |VectorArray| containing a basis of the reduced space onto which to project.
    product
        Inner product |Operator| w.r.t. which to compute the Riesz representatives
        in case `rhs_is_functional` is `True`. When `product` is `None`, no Riesz
        representatives are computed
    extends
        Set by :meth:`~pymor.algorithms.greedy.greedy` to the result of the
        last reduction in case the basis extension was `hierarchic` (used to prevent
        re-computation of `residual_range` basis vectors already obtained from previous
        reductions).

    Returns
    -------
    projected_residual
        See above.
    reconstructor
        See above.
    reduction_data
        Additional data produced by the reduction process (compare the `extends` parameter).
    """
    assert rhs is None \
        or rhs_is_functional and (rhs.range == NumpyVectorSpace(1) and rhs.source == operator.range and rhs.linear) \
        or not rhs_is_functional and (rhs.source == NumpyVectorSpace(1) and rhs.range == operator.range and rhs.linear)
    assert RB is None or RB in operator.source
    assert product is None or product.source == product.range == operator.range
    assert extends is None or len(extends) == 3

    logger = getLogger('pymor.reductors.residual.reduce_residual')

    if RB is None:
        RB = operator.source.empty()

    if extends and isinstance(extends[0], NonProjectedResidualOperator):
        extends = None
    if extends:
        residual_range = extends[1].RB
        residual_range_dims = list(extends[2]['residual_range_dims'])
    else:
        residual_range = operator.range.empty()
        residual_range_dims = []

    with logger.block('Estimating residual range ...'):
        try:
            residual_range, residual_range_dims = \
                estimate_image_hierarchical([operator], [rhs], RB, (residual_range, residual_range_dims),
                                            orthonormalize=True, product=product,
                                            riesz_representatives=rhs_is_functional)
        except ImageCollectionError as e:
            logger.warn('Cannot compute range of {}. Evaluation will be slow.'.format(e.op))
            operator = operator.projected(None, RB)
            return (NonProjectedResidualOperator(operator, rhs, rhs_is_functional, product),
                    NonProjectedReconstructor(product),
                    {})

    with logger.block('Projecting residual operator ...'):
        if rhs_is_functional:
            operator = operator.projected(residual_range, RB, product=None)  # the product cancels out.
            rhs = rhs.projected(None, residual_range, product=None)
        else:
            operator = operator.projected(residual_range, RB, product=product)
            rhs = rhs.projected(residual_range, None, product=product)

    return (ResidualOperator(operator, rhs, rhs_is_functional),
            GenericRBReconstructor(residual_range),
            {'residual_range_dims': residual_range_dims})

class ResidualOperator(OperatorBase):
    """Returned by :func:`reduce_residual`."""

    def __init__(self, operator, rhs, rhs_is_functional=True, name=None):
        self.source = operator.source
        self.range = operator.range
        self.linear = operator.linear
        self.operator = operator
        self.rhs = rhs
        self.rhs_vector = rhs.as_vector() if rhs and not rhs.parametric else None
        self.rhs_is_functional = rhs_is_functional
        self.name = name

    def apply(self, U, ind=None, mu=None):
        V = self.operator.apply(U, ind=ind, mu=mu)
        if self.rhs:
            F = self.rhs_vector or self.rhs.as_vector(mu)
            if len(V) > 1:
                V.axpy(-1., F, x_ind=[0]*len(V))
            else:
                V.axpy(-1., F)
        return V

    def projected_to_subbasis(self, dim_range=None, dim_source=None, name=None):
        if self.rhs_is_functional:
            rhs = self.rhs.projected_to_subbasis(None, dim_range)
        else:
            rhs = self.rhs.projected_to_subbasis(dim_range, None)
        return ResidualOperator(self.operator.projected_to_subbasis(dim_range, dim_source), rhs,
                                self.rhs_is_functional, name=name)

class NonProjectedResidualOperator(ResidualOperator):
    """Returned by :func:`reduce_residual`.

    Not to be used directly.
    """

    def __init__(self, operator, rhs, rhs_is_functional, product):
        super().__init__(operator, rhs, rhs_is_functional)
        self.product = product

    def apply(self, U, ind=None, mu=None):
        R = super().apply(U, ind=ind, mu=mu)
        if self.product:
            if self.rhs_is_functional:
                R_riesz = self.product.apply_inverse(R)
                return R_riesz * (np.sqrt(R_riesz.dot(R)) / R_riesz.l2_norm())[0]
            else:
                return R * (np.sqrt(self.product.pairwise_apply2(R, R)) / R.l2_norm())[0]
        else:
            return R

    def projected_to_subbasis(self, dim_range=None, dim_source=None, name=None):
        return self.with_(operator=self.operator.projected_to_subbasis(None, dim_source))

class NonProjectedReconstructor(ImmutableInterface):
    """Returned by :func:`reduce_residual`.

    Not to be used directly.
    """

    def __init__(self, product):
        self.norm = induced_norm(product) if product else None

    def reconstruct(self, U):
        if self.norm:
            return U * (U.l2_norm() / self.norm(U))[0]
        else:
            return U

def reduce_implicit_euler_residual(operator, mass, dt, functional=None, RB=None, product=None, extends=None):
    """Reduced basis residual reductor with mass operator for implicit Euler timestepping.

    Given an operator, mass and a functional, the concatenation of residual operator
    with the Riesz isomorphism is given by::

        riesz_residual.apply(U, U_old, mu)
            == product.apply_inverse(operator.apply(U, mu) + 1/dt*mass.apply(U, mu) - 1/dt*mass.apply(U_old, mu)
               - functional.as_vector(mu))

    This reductor determines a low-dimensional subspace of the image of a reduced
    basis space under `riesz_residual` using :func:`~pymor.algorithms.image.estimate_image_hierarchical`,
    computes an orthonormal basis `residual_range` of this range space and then
    returns the Petrov-Galerkin projection ::

        projected_riesz_residual
            === riesz_residual.projected(range_basis=residual_range, source_basis=RB)

    of the `riesz_residual` operator. Given reduced basis coefficient vectors `u` and `u_old`,
    the dual norm of the residual can then be computed as ::

        projected_riesz_residual.apply(u, u_old, mu).l2_norm()

    Moreover, a `reconstructor` is provided such that ::

        reconstructor.reconstruct(projected_riesz_residual.apply(u, u_old, mu))
            == riesz_residual.apply(RB.lincomb(u), RB.lincomb(u_old), mu)

    Parameters
    ----------
    operator
        See definition of `riesz_residual`.
    mass
        The mass operator. See definition of `riesz_residual`.
    dt
        The time step size. See definition of `riesz_residual`.
    functional
        See definition of `riesz_residual`. If `None`, zero right-hand side is assumed.
    RB
        |VectorArray| containing a basis of the reduced space onto which to project.
    product
        Inner product |Operator| w.r.t. which to compute the Riesz representatives.
    extends
        Set by :meth:`~pymor.algorithms.greedy.greedy` to the result of the
        last reduction in case the basis extension was `hierarchic` (used to prevent
        re-computation of `residual_range` basis vectors already obtained from previous
        reductions).

    Returns
    -------
    projected_riesz_residual
        See above.
    reconstructor
        See above.
    reduction_data
        Additional data produced by the reduction process (compare the `extends` parameter).
    """
    assert functional is None \
        or functional.range == NumpyVectorSpace(1) and functional.source == operator.source and functional.linear
    assert RB is None or RB in operator.source
    assert product is None or product.source == product.range == operator.range
    assert extends is None or len(extends) == 3

    logger = getLogger('pymor.reductors.residual.reduce_implicit_euler_residual')

    if RB is None:
        RB = operator.source.empty()

    if extends and isinstance(extends[0], NonProjectedImplicitEulerResidualOperator):
        extends = None
    if extends:
        residual_range = extends[1].RB
        residual_range_dims = list(extends[2]['residual_range_dims'])
    else:
        residual_range = operator.range.empty()
        residual_range_dims = []

    with logger.block('Estimating residual range ...'):
        try:
            residual_range, residual_range_dims = \
                estimate_image_hierarchical([operator, mass], [functional], RB, (residual_range, residual_range_dims),
                                            orthonormalize=True, product=product, riesz_representatives=True)
        except ImageCollectionError as e:
            logger.warn('Cannot compute range of {}. Evaluation will be slow.'.format(e.op))
            operator = operator.projected(None, RB)
            mass = mass.projected(None, RB)
            return (NonProjectedImplicitEulerResidualOperator(operator, mass, functional, dt, product),
                    NonProjectedReconstructor(product),
                    {})

    with logger.block('Projecting residual operator ...'):
        operator = operator.projected(residual_range, RB, product=None)  # the product always cancels out.
        mass = mass.projected(residual_range, RB, product=None)
        functional = functional.projected(None, residual_range, product=None)

    return (ImplicitEulerResidualOperator(operator, mass, functional, dt),
            GenericRBReconstructor(residual_range),
            {'residual_range_dims': residual_range_dims})

class ImplicitEulerResidualOperator(OperatorBase):
    """Returned by :func:`reduce_implicit_euler_residual`."""

    def __init__(self, operator, mass, functional, dt, name=None):
        self.source = operator.source
        self.range = operator.range
        self.linear = operator.linear
        self.operator = operator
        self.mass = mass
        self.functional = functional
        self.functional_vector = functional.as_vector() if functional and not functional.parametric else None
        self.dt = dt
        self.name = name

    def apply(self, U, U_old, ind, ind_old=None, mu=None):
        V = self.operator.apply(U, ind=ind, mu=mu)
        V.axpy(1./self.dt, self.mass.apply(U, ind=ind, mu=mu))
        V.axpy(-1./self.dt, self.mass.apply(U_old, ind=ind_old, mu=mu))
        if self.functional:
            F = self.functional_vector or self.functional.as_vector(mu)
            if len(V) > 1:
                V.axpy(-1., F, x_ind=[0]*len(V))
            else:
                V.axpy(-1., F)
        return V

    def projected_to_subbasis(self, dim_range=None, dim_source=None, name=None):
        return ImplicitEulerResidualOperator(self.operator.projected_to_subbasis(dim_range, dim_source),
                                             self.mass.projected_to_subbasis(dim_range, dim_source),
                                             self.functional.projected_to_subbasis(None, dim_range),
                                             self.dt,
                                             name=name)

class NonProjectedImplicitEulerResidualOperator(ImplicitEulerResidualOperator):
    """Returned by :func:`reduce_implicit_euler_residual`.

    Not to be used directly.
    """

    def __init__(self, operator, mass, functional, dt, product):
        super().__init__(operator, mass, functional, dt)
        self.product = product

    def apply(self, U, U_old, ind=None, ind_old=None, mu=None):
        R = super().apply(U, U_old, ind=ind, ind_old=ind_old, mu=mu)
        if self.product:
            R_riesz = self.product.apply_inverse(R)
            return R_riesz * (np.sqrt(R_riesz.dot(R)) / R_riesz.l2_norm())[0]
        else:
            return R

    def projected_to_subbasis(self, dim_range=None, dim_source=None, name=None):
        return self.with_(operator=self.operator.projected_to_subbasis(None, dim_source),
                          mass=self.mass.projected_to_subbasis(None, dim_source))

max_stars_repo_name: JuliaBru/pymor

max_issues_repo_name: JuliaBru/pymor

max_forks_repo_name: JuliaBru/pymor


text:
import numpy as np
import tensorflow as tf
from tensorflow.keras import Model
import matplotlib.pyplot as plt

# preprocessing
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import MinMaxScaler
from sklearn.preprocessing import StandardScaler

# Some metrics
# Coefficient of determination
def coeff_determination(y_pred, y_true): #Order of function inputs is important here
    SS_res =  np.sum(np.square( y_true-y_pred ))
    SS_tot = np.sum(np.square( y_true - np.mean(y_true) ) )
    return ( 1 - SS_res/(SS_tot + 2.22044604925e-16) )

#Build the model which does basic map of inputs to coefficients
class standard_lstm(Model):
    def __init__(self,data):
        super(standard_lstm, self).__init__()

        # Set up the data for the LSTM
        self.data_tsteps = np.shape(data)[0]
        self.state_len = np.shape(data)[1]

        self.preproc_pipeline = Pipeline([('stdscaler', StandardScaler()),('minmax', MinMaxScaler(feature_range=(-1, 1)))])
        self.data = self.preproc_pipeline.fit_transform(data)

        # Need to make minibatches
        self.seq_num = 5
        self.total_size = np.shape(data)[0]-int(self.seq_num) # Limit of sampling

        input_seq = np.zeros(shape=(self.total_size,self.seq_num,self.state_len))  #[samples,n_inputs,state_len]
        output_seq = np.zeros(shape=(self.total_size,self.state_len)) #[samples,n_outputs,state_len]

        snum = 0
        for t in range(0,self.total_size):
            input_seq[snum,:,:] = self.data[None,t:t+self.seq_num,:]
            output_seq[snum,:] = self.data[None,t+self.seq_num,:]
            snum = snum + 1

        # Shuffle dataset
        idx = np.arange(snum)
        np.random.shuffle(idx)
        input_seq = input_seq[idx]
        output_seq = output_seq[idx]

        # Split into train and test
        self.input_seq_test = input_seq[int(0.9*snum):]
        self.output_seq_test = output_seq[int(0.9*snum):]
        input_seq = input_seq[:int(0.9*snum)]
        output_seq = output_seq[:int(0.9*snum)]

        # Split into train and valid
        self.ntrain = int(0.8*np.shape(input_seq)[0])
        self.nvalid = np.shape(input_seq)[0] - self.ntrain

        self.input_seq_train = input_seq[:self.ntrain]
        self.output_seq_train = output_seq[:self.ntrain]

        self.input_seq_valid = input_seq[self.ntrain:]
        self.output_seq_valid = output_seq[self.ntrain:]

        # Define architecture
        xavier=tf.keras.initializers.GlorotUniform()
        self.l1=tf.keras.layers.LSTM(50,return_sequences=True,input_shape=(self.seq_num,self.state_len))
        self.l2=tf.keras.layers.LSTM(50,return_sequences=False)
        self.out = tf.keras.layers.Dense(self.state_len)
        self.train_op = tf.keras.optimizers.Adam(learning_rate=0.001)

    # Running the model
    def call(self,X):
        h1 = self.l1(X)
        h2 = self.l2(h1)
        out = self.out(h2)
        return out

    # Regular MSE
    def get_loss(self,X,Y):
        op=self.call(X)
        return tf.reduce_mean(tf.math.square(op-Y))

    # get gradients - regular
    def get_grad(self,X,Y):
        with tf.GradientTape() as tape:
            tape.watch(self.trainable_variables)
            L = self.get_loss(X,Y)
            g = tape.gradient(L, self.trainable_variables)
        return g

    # perform gradient descent - regular
    def network_learn(self,X,Y):
        g = self.get_grad(X,Y)
        self.train_op.apply_gradients(zip(g, self.trainable_variables))

    # Train the model
    def train_model(self):
        plot_iter = 0
        stop_iter = 0
        patience = 10
        best_valid_loss = 999999.0 # Some large number

        self.num_batches = 20
        self.train_batch_size = int(self.ntrain/self.num_batches)
        self.valid_batch_size = int((self.nvalid)/self.num_batches)

        for i in range(10):
            # Training loss
            print('Training iteration:',i)

            for batch in range(self.num_batches):
                input_batch = self.input_seq_train[batch*self.train_batch_size:(batch+1)*self.train_batch_size]
                output_batch = self.output_seq_train[batch*self.train_batch_size:(batch+1)*self.train_batch_size]
                self.network_learn(input_batch,output_batch)

            # Validation loss
            valid_loss = 0.0
            valid_r2 = 0.0

            for batch in range(self.num_batches):
                input_batch = self.input_seq_valid[batch*self.valid_batch_size:(batch+1)*self.valid_batch_size]
                output_batch = self.output_seq_valid[batch*self.valid_batch_size:(batch+1)*self.valid_batch_size]

                valid_loss = valid_loss + self.get_loss(input_batch,output_batch).numpy()
                predictions = self.call(self.input_seq_valid)
                valid_r2 = valid_r2 + coeff_determination(predictions,self.output_seq_valid)

            valid_r2 = valid_r2/(batch+1)

            # Check early stopping criteria
            if valid_loss < best_valid_loss:

                print('Improved validation loss from:',best_valid_loss,' to:', valid_loss)
                print('Validation R2:',valid_r2)

                best_valid_loss = valid_loss

                self.save_weights('./checkpoints/my_checkpoint')

                stop_iter = 0
            else:
                print('Validation loss (no improvement):',valid_loss)
                print('Validation R2:',valid_r2)
                stop_iter = stop_iter + 1

            if stop_iter == patience:
                break

        # Check accuracy on test
        predictions = self.call(self.input_seq_test)
        print('Test loss:',self.get_loss(self.input_seq_test,self.output_seq_test).numpy())
        r2 = coeff_determination(predictions,self.output_seq_test)
        print('Test R2:',r2)
        r2_iter = 0

    # Load weights
    def restore_model(self):
        self.load_weights('./checkpoints/my_checkpoint') # Load pretrained model
        print('Model restored successfully!')

    # Do some testing
    def model_inference(self,test_data):
        # Restore from checkpoint
        self.restore_model()

        # Scale testing data
        test_data = self.preproc_pipeline.fit_transform(test_data)

        # Test sizes
        test_total_size = np.shape(test_data)[0]-int(self.seq_num) # Limit of sampling

        # Input placeholder setup
        rec_input_seq = test_data[:self.seq_num,:].reshape(1,self.seq_num,self.state_len)

        # True outputs
        rec_output_seq = np.zeros(shape=(test_total_size,self.state_len)) #[samples,n_outputs,state_len]
        snum = 0
        for t in range(test_total_size):
            rec_output_seq[snum,:] = test_data[None,t+self.seq_num,:]
            snum = snum + 1

        # Recursive predict
        print('Making predictions on testing data')
        rec_pred = np.copy(rec_output_seq)
        for t in range(test_total_size):
            rec_pred[t] = self.call(rec_input_seq).numpy()[0]
            rec_input_seq[0,0:-1,:] = rec_input_seq[0,1:]
            rec_input_seq[0,-1,:] = rec_pred[t]

        # Rescale
        rec_pred = self.preproc_pipeline.inverse_transform(rec_pred)
        rec_output_seq = self.preproc_pipeline.inverse_transform(rec_output_seq)

        # plot
        for i in range(self.state_len):
            plt.figure()
            plt.title('Mode '+str(i))
            plt.plot(rec_pred[:,i],label='Predicted')
            plt.plot(rec_output_seq[:,i],label='True')
            plt.legend()
            plt.savefig('Mode_'+str(i)+'_prediction.png')
            plt.close()

        return rec_output_seq, rec_pred

if __name__ == '__main__':
    print('Architecture file')

max_stars_repo_name: Romit-Maulik/Tutorials-Demos-Practice

max_issues_repo_name: omersan/Practice

max_forks_repo_name: omersan/Practice


text:
import cv2
import numpy as np

def apply_hist_mask(frame, hist):
    hsv = cv2.cvtColor(frame, cv2.COLOR_BGR2HSV)
    dst = cv2.calcBackProject([hsv], [0,1], hist, [0,180,0,256], 1)

    disc = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (11,11))
    cv2.filter2D(dst, -1, disc, dst)

    ret, thresh = cv2.threshold(dst, 100, 255, 0)
    thresh = cv2.merge((thresh,thresh, thresh))

    cv2.GaussianBlur(dst, (3,3), 0, dst)

    res = cv2.bitwise_and(frame, thresh)
    cv2.imshow("images", res)
    cv2.moveWindow("images", 0, 600)
    return res

def contours(frame):
    gray = cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY)
    ret,thresh = cv2.threshold(gray, 0, 255, 0)
    frame, contours, hierarchy = cv2.findContours(thresh, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
    return contours

def max_contour(contours):
    max_i = 0
    max_area = 0

    for i in range(len(contours)):
        cnt = contours[i]
        area = cv2.contourArea(cnt)
        if area > max_area:
            max_area = area
            max_i = i

    contour = contours[max_i]
    return contour

def hull(contour):
    hull = cv2.convexHull(contour)
    return hull

def defects(contour):
    hull = cv2.convexHull(contour, returnPoints=False)
    if hull is not None and len(hull > 3) and len(contour) > 3:
        defects = cv2.convexityDefects(contour, hull)
        return defects
    else:
        return None

def centroid(contour):
    moments = cv2.moments(contour)
    if moments['m00'] != 0:
        cx = int(moments['m10']/moments['m00'])
        cy = int(moments['m01']/moments['m00'])
        return (cx,cy)
    else:
        return None

def contour_interior(frame, contour):
    rect = cv2.minAreaRect(contour)
    box = cv2.cv.BoxPoints(rect)
    box = np.int0(box)

    rows,cols,_ = frame.shape
    mask = np.zeros((rows,cols), dtype=np.float)
    for i in xrange(rows):
        for j in xrange(cols):
            mask.itemset((i,j), cv2.pointPolygonTest(box, (j,i), False))

    mask = np.int0(mask)
    mask[mask < 0] = 0
    mask[mask > 0] = 255
    mask = np.array(mask, dtype=frame.dtype)
    mask = cv2.merge((mask, mask, mask))

    contour_interior = cv2.bitwise_and(frame, mask)
    return contour_interior

def gray_threshold(frame, threshold_value):
    gray = cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY)
    ret, thresh = cv2.threshold(gray, threshold_value, 255, 0)
    return thresh

def farthest_point(defects, contour, centroid):
    s = defects[:,0][:,0]
    cx, cy = centroid

    x = np.array(contour[s][:,0][:,0], dtype=np.float)
    y = np.array(contour[s][:,0][:,1], dtype=np.float)

    xp = cv2.pow(cv2.subtract(x, cx), 2)
    yp = cv2.pow(cv2.subtract(y, cy), 2)
    dist = cv2.sqrt(cv2.add(xp, yp))

    dist_max_i = np.argmax(dist)

    if dist_max_i < len(s):
        farthest_defect = s[dist_max_i]
        farthest_point = tuple(contour[farthest_defect][0])
        return farthest_point
    else:
        return None


text:
import numpy as np

from time import time
from qiskit.circuit import QuantumCircuit, ParameterVector
from qiskit.circuit.library import RealAmplitudes

from surfer.qfi import StochasticApproximation, LinearCombination, ReverseQFI

circuit = RealAmplitudes(5, reps=15)
values = np.ones(circuit.num_parameters)
print(circuit.num_parameters)
# x = ParameterVector('x', 2)
# circuit = QuantumCircuit(1)
# circuit.ry(x[0], 0)
# circuit.ry(x[1], 0)
# values = np.array([0.5, 0])

# reference implementation
# lc = LinearCombination()
# start = time()
# reference = lc.compute(circuit, values)

# print(f"Reference: {time() - start}")
# print(reference)

# stochastic approximation
# sa = StochasticApproximation(samples=1000, perturbation=0.1)
# approximated = sa.compute(circuit, values)

# print('Approximated:')
# print(approximated)

rev = ReverseQFI()
start = time()
revd = rev.compute(circuit, values)

print(f"Reverse mode: {time() - start}")
# print(np.mean(np.abs(revd - reference)))

max_stars_repo_name: Cryoris/surfer

max_issues_repo_name: Cryoris/surfer

max_forks_repo_name: Cryoris/surfer


text:
# python src/chapter10/chapter10note.py
# python3 src/chapter10/chapter10note.py
'''
Class Chapter11_1

'''
import sys as _sys
import math as _math
import random as _random
import time as _time
from random import randint as _randint
from copy import copy as _copy, deepcopy as _deepcopy
from numpy import arange as _arange
import numpy as np

class Chapter11_1:
    '''
    chpater11.1 note and function
    '''
    def note(self):
        '''
        Summary
        ====
        Print chapter11.1 note

        Example
        ====
        ```python
        Chapter11_1().note()
        ```
        '''
        print('chapter11.1 note as follow')
        print('第11章 散列表')
        print('在很多应用中，都要用到一种动态集合结构，它仅支持INSERT,SEARCH的DELETE字典操作')
        print('实现字典的一种有效数据结构为散列表(HashTable)')
        print('在最坏情况下，在散列表中，查找一个元素的时间在与链表中查找一个元素的时间相同')
        print('在最坏情况下都是Θ(n)')
        print('但在实践中，散列技术的效率是很高的。在一些合理的假设下，在散列表中查找一个元素的期望时间为O(1)')
        print('散列表是普通数组概念的推广，因为可以对数组进行直接寻址，故可以在O(1)时间内访问数组的任意元素')
        print('如果存储空间允许，可以提供一个数组，为每个可能的关键字保留一个位置，就可以应用直接寻址技术')
        print('当实际存储的关键字数比可能的关键字总数较小时，这是采用散列表就会较直接数组寻址更为有效')
        print('在散列表中，不是直接把关键字用作数组下标，而是根据关键字计算出下标。')
        print('11.2着重介绍解决碰撞的链接技术。')
        print('所谓碰撞，就是指多个关键字映射到同一个数组下标位置')
        print('11.3介绍如何利用散列函数，根据关键字计算出数组的下标。')
        print('11.4介绍开放寻址法，它是处理碰撞的另一种方法。散列是一种极其有效和实用的技术，基本的字典操作只需要O(1)的平均时间')
        print('11.5解释当待排序的关键字集合是静态的，\"完全散列\"如何能够在O(1)最坏情况时间内支持关键字查找')
        print('11.1 直接寻址表')
        print('当关键字的全域U比较小时，直接寻址是一种简单而有效的技术。假设某应用要用到一个动态集合')
        print('其中每个元素都有一个取自全域U的关键字，此处m是一个不很大的数。另外假设没有两个元素具有相同的关键字')
        print('为表示动态集合，用一个数组(或直接寻址表)T[0...m-1],其中每个位置(或称槽)对应全域U中的一个关键字')
        print('对于某些应用，动态集合中的元素可以放在直接寻址表中。亦即不把每个元素的关键字及其卫星数据都放在直接寻址表外部的一个对象中')
        print('但是，如果不存储关键字，就必须有某种办法来确定某个槽是否为空')
        print('练习11.1-1: 考虑一个由长度为m的直接寻址表T表示的动态集合S。给出一个查找S的最大元素的算法过程')
        print(' 所给的过程在最坏情况下的运行时间是O(m)')
        print('练习11.1-2: 位向量(bit vector)是一种仅包含0和1的数组。长度为m的位向量所占空间要比包含m个指针的数组少得多')
        print(' 请说明如何用一个位向量来表示一个包含不同元素的动态集合。字典操作的运行时间应该是O(1)')
        print('练习11.1-3: 说明如何实现一个直接寻址表，使各元素的关键字不必都相同，且各元素可以有卫星数据。')
        print('练习11.1-4: 希望通过利用一个非常大的数组上直接寻址的方式来实现字典')
        print(' 开始时，该数组中可能包含废料，但要对整个数组进行初始化是不实际的，因为该组的规模太大')
        print(' 请给出在大数组上实现直接寻址字典的方案。每个存储的对象占用O(1)空间')
        print(' 操作SEARCH,INSERT和DELETE的时间为O(1),对数据结构初始化的时间O(1)')
        print(' 可以利用另外一个栈，其大小等于实际存储在字典中的关键字数目，以帮助确定大型数组中某个给定的项是否是有效的')
        # python src/chapter11/chapter11note.py
        # python3 src/chapter11/chapter11note.py

class Chapter11_2:
    '''
    chpater11.2 note and function
    '''
    def note(self):
        '''
        Summary
        ====
        Print chapter11.2 note

        Example
        ====
        ```python
        Chapter11_2().note()
        ```
        '''
        print('chapter11.2 note as follow')
        print('11.2 散列表')
        print('直接寻址技术存在着一个明显的问题：如果域U很大，',
            '在一台典型计算机的可用内存容量限制下，要在机器中存储大小为U的一张表T就有点不实际甚至是不可能的了')
        print('实际要存储的关键字集合K相对于U来说可能很小，因而分配给T的大部分空间都要浪费掉')
        print('当存储在字典中的关键字集合K比所有可能的关键字域U要小的多时，散列表需要的存储空间要比直接寻址少很多')
        print('特别地，在保持仅需O(1)时间即可在散列表中查找一个元素的好处情况下，存储要求可以降至Θ(|K|)')
        print('在直接寻址方式下，具有关键字k的元素被存放在槽k中。在散列方式下，该元素处于h(k)中')
        print('亦即，利用散列函数h,根据关键字k计算出槽的位置。函数h将关键字域U映射懂啊散列表T[0..m-1]的槽位上：')
        print('这时，可以说一个具有关键字k到元素是被散列在槽h(k)上，或说h(k)是关键字k的散列值')
        print('两个关键字可能映射到同一个槽上。将这种情形称为发生了碰撞')
        print('当然，最理想的解决方法是完全避免碰撞')
        print('可以考虑选用合适的散列函数h。在选择时有一个主导思想，就是使h尽可能地\"随机\",从而避免或至少最小化碰撞')
        print('当然，一个散列函数h必须是确定的，即某一给定的输入k应始终产生相同的结果h(k),')
        print('通过链接法解决碰撞')
        print(' 在链接法中，把散列到同一槽中的所有元素都放在同一个链表中，槽j中有一个指针，它指向由所有散列到j的元素构成的链表的头')
        print(' 插入操作的最坏情况运行时间为O(1).插入过程要快一些，因为假设要插入的元素x没有出现在表中；如果需要，在插入前执行搜索，可以检查这个假设(付出额外代价)')
        print('CHAINED-HASH-INSERT(T, x)')
        print(' insert x at the head of list T[h(key[x])]')
        print('CHAINED-HASH-SEARCH(T, x)')
        print(' search for an element with k in list T[h(k)]')
        print('CHAINED-HASH-DELETE(T, x)')
        print(' delete x from the list T[h(key[x])]')
        print('对用链接法散列的分析')
        print(' 采用链接法后散列的性能怎样呢？特别地，要查找一个具有给定关键字的原宿需要多长时间呢？')
        print(' 给定一个能存放n个元素的，具有m个槽位的散列表T，定义T的装载因子a为n/m,即一个链表中平均存储的元素数')
        print(' 分析以a来表达，a可以小于、等于或者大于1')
        print('用链接法散列的最坏情况性能很差；所有的n个关键字都散列到同一个槽中，从而产生出一个长度为n的链表')
        print('最坏情况下查找的时间为Θ(n),再加上计算散列函数的时间，这么一来就和用一个链表来来链接所有的元素差不多了。显然，')
        print('散列方法的平均性态依赖于所选取的散列函数h在一般情况下，将所有的关键字分布在m个槽位上的均匀程度')
        print('先假定任何元素散列到m个槽中每一个的可能性是相同的，且与其他元素已被散列到什么位置上是独立无关的')
        print('称这个假设为简单一致散列')
        print('假定可以在O(1)时间内计算出散列值h(k),从而查找具有关键字为k的元素的时间线性地依赖于表T[h(k)]的长度为n')
        print('先不考虑计算散列函数和寻址槽h(k)的O(1)时间')
        print('定理11.1 对一个用链接技术来解决碰撞的散列表，在简单一致散列的假设下，一次不成功查找期望时间为Θ(1+a)')
        print('定理11.2 在简单一致散列的假设下，对于用链接技术解决碰撞的散列表，平均情况下一次成功的查找需要Θ(1+a)时间')
        print('练习11.2-1: 假设用一个散列函数h，将n个不同的关键字散列到一个长度为m的数组T中。')
        print(' 假定采用的是简单一致散列法，那么期望的碰撞数是多少？')
        print('练习11.2-2: 对于一个利用链接法解决碰撞的散列表，说明将关键字5,28,19,15,20,33,12,17,10')
        print(' 设该表中有9个槽位，并设散列函数为h(k)=k mod 9')
        print('练习11.2-3: 如果将链接模式改动一下，使得每个链表都能保持已排序顺序，散列的性能就可以有很大的提高。')
        print(' 这样的改动对成功查找、不成功查找、插入和删除操作的运行时间有什么影响')
        print('练习11.2-4: 在散列表内部，如何通过将所有未占用的槽位链接成一个自由链表，来分配和去分配元素的存储空间')
        print(' 假定一个槽位可以存储一个标志、一个元素加上一个或两个指针')
        print(' 所有的字典和自由链表操作应具有O(1)的期望运行时间')
        print('练习11.2-5: 有一个U的大小为n的子集，它包含了均散列到同一个槽位中的关键字，这样对于带链接的散列表，最坏情况下查找时间为Θ(n)')
        # python src/chapter11/chapter11note.py
        # python3 src/chapter11/chapter11note.py

class Chapter11_3:
    '''
    chpater11.3 note and function
    '''
    def note(self):
        '''
        Summary
        ====
        Print chapter11.3 note

        Example
        ====
        ```python
        Chapter11_3().note()
        ```
        '''
        print('chapter11.3 note as follow')
        print('11.3 散列函数')
        print('好的散列函数的特点')
        print('一个好的散列函数应近似地满足简单一致散列的假设：每个关键字都等可能地散列到m个槽位的任何一个之中去')
        print('并与其它的关键字已被散列到哪一个槽位中无关')
        print('不幸的是：一般情况下不太可能检查这一条件是否成立，因为人们很少可能知道关键字所符合的概率分布，而各关键字可能并不是完全互相独立的')
        print('有时也能知道关键字的概率分布。例如：已知各关键字都是随机的实数k，独立地、一致地分布于范围[0,1)')
        print('在实践中，常常可以运用启发式技术来构造好的散列函数')
        print('例如，在一个编译器的符号表中，关键字都是字符串，表示程序中的标识符')
        print('同一个程序中，经常会出现一些很相近的符号，如pt和pts。')
        print('一个好的散列函数应能最小化将这些相近符号散列到同一个槽中的可能性')
        print('\"除法散列\"用一个特定的质数来除所给的关键字，所得的余数即为该关键字的散列值')
        print('假定所选择的质数与关键字分布中的任何模式都是无关的，这种方法常常可以给出很好的结果')
        print('散列函数的某些应用可能会要求比简单一致散列更强的性质，例如可能希望某些很近似的关键字具有截然不同的散列值')
        print('将关键字解释为自然数')
        print(' 如果所给关键字不是自然数，则必须有一种方法来将它们解释为自然数')
        print(' 标识符pt可以被解释为十进制整数对(112,116),pt即为(112*128)+116=14452')
        print('11.3.1 除法散列法')
        print(' 通过取k除以m的余数，来将关键字k映射到m个槽的摸一个中去，亦即散列函数为h(k) = k mod m')
        print(' 例如，如果散列表的大小为m=12,所给关键字为k=100,则h(k)=4。这种方法只要一次除法操作，所以比较快')
        print(' 应用除法散列时，要注意m的选择，m不应是2的幂；可以选作m的值常常是与2的整数幂不太接近的质数')
        print(' 例如，假设我们要分配一张散列表，并用链接法解决碰撞，表中大约要存放n=2000个字符串，每个字符有8位')
        print(' 一次不成功的查找大约要检查3个元素，但我们并不在意，故分配散列表的大小为m=701.')
        print(' 之所以选择701这个数，是因为它是个接近a=2000/3，但又不接近2的任何幂次的质数。把每个关键字k视为一个整数')
        print(' 则有散列函数h(k) = k mod 701')
        print('11.3.2 乘法散列法')
        print(' 构造散列函数的乘法方法包含两个步骤。第一步，用关键字k乘上常数A(0<A<1),',
            '并抽出kA的小数部分。然后，用m乘以这个值，再取结果的底(floor)。总之，散列函数为h(k)=[m(kA mod 1)]')
        print(' 其中kA mod 1 即kA的小数部分，亦即kA-[kA]')
        print(' 乘法方法的一个优点是对m的选择没有特别的要求，一般选择它为2的某个幂次m=2^p')
        print(' 虽然这个方法对任何的A值都适用，但对某些值效果更好。Knuth认为最佳的选择与待散列的数据的特征有关A=(sqrt(5)-1)/2=0.6180339887...就是一个比较理想的值')
        print(' 例子：假设有k=123456,p=14,m=2^14=16384,w=32,根据Knuth的建议')
        print(' 取A为形如s/2^32的分数，它与(sqrt(5)-1)/2最为接近，于是A=2654435769/2^32')
        print(' k*s=32770622297664=(76300*2^32)+17612864,从而有r1=76300,r0=17612864,r0的14个最高有效位产生了散列值h(k)=67')
        print('11.3.3 全域散列')
        print('如果让某个与你作对的人来选择要散列的关键字，那么他会选择全部散列到同一槽中的n个关键字，使得平均检索值为Θ(n)')
        print('任何一个特定的散列函数都可能出现这种最坏情况性态：唯一有效的改进方法是随机地选择散列函数，使之独立要存储的关键字。')
        print('这种方法称作全域散列(universal hashing),不管对手选择了怎样的关键字，其平均性态都很好')
        print('全域散列的基本思想是在执行开始时，就从一族仔细设计的函数中，随机地选择一个座位散列函数')
        print('就像在快速排序中一样，随机化保证了没有哪一种输入会导致最坏情况性态。')
        print('同时，随机化使得即使对同一个输入，算法在每一次执行时的性态也都不一样')
        print('这样就可以确保对于任何输入，算法都具有较好的平均情况性态')
        print('设H为有限的一组散列函数，它将给定的关键字域U映射到{0,1,..,m-1}中。这样的一个函数组称为是全域的')
        print('定理11.3 如果h选自一组全域的散列函数，并用于将n个关键字散列到一个大小为m的、用链接法解决碰撞的表T中')
        print('对于每一个关键字k，定义一个随机变量Yk,它等于非k的、与k散列到同一槽位中的其他关键字的数目')
        print('推论11.4 对于一个具有m个槽位的表，利用全域散列和链接法解决碰撞，需要Θ(n)的期望时间来处理任何包含了n')
        print(' 个INSERT,SEARCH和DELETE操作的操作序列，该序列中包含了O(m)个INSERT操作')
        print('证明：由于插入操作的数目为O(m),有n=O(m),从而a=O(1),')
        print('')
        print('')
        # python src/chapter11/chapter11note.py
        # python3 src/chapter11/chapter11note.py

class Chapter11_4:
    '''
    chpater11.4 note and function
    '''
    def note(self):
        '''
        Summary
        ====
        Print chapter11.4 note

        Example
        ====
        ```python
        Chapter11_4().note()
        ```
        '''
        print('chapter11.4 note as follow')
        # python src/chapter11/chapter11note.py
        # python3 src/chapter11/chapter11note.py

class Chapter11_5:
    '''
    chpater11.5 note and function
    '''
    def note(self):
        '''
        Summary
        ====
        Print chapter11.5 note

        Example
        ====
        ```python
        Chapter11_5().note()
        ```
        '''
        print('chapter11.5 note as follow')
        # python src/chapter11/chapter11note.py
        # python3 src/chapter11/chapter11note.py

chapter11_1 = Chapter11_1()
chapter11_2 = Chapter11_2()
chapter11_3 = Chapter11_3()
chapter11_4 = Chapter11_4()
chapter11_5 = Chapter11_5()

def printchapter11note():
    '''
    print chapter11 note.
    '''
    print('Run main : single chapter eleven!')
    chapter11_1.note()
    chapter11_2.note()
    chapter11_3.note()
    chapter11_4.note()
    chapter11_5.note()

# python src/chapter10/chapter10note.py
# python3 src/chapter10/chapter10note.py
if __name__ == '__main__':
    printchapter11note()
else:
    pass

reason: import numpy,from numpy


text:
#!/usr/bin/env python3
"""Imports nematus model file and convert it into a neural monkey experiment
given a neural monkey configuration file.
"""
from typing import Dict, Tuple, List
import argparse
import json
import os
import numpy as np
import tensorflow as tf

from neuralmonkey.config.parsing import parse_file
from neuralmonkey.config.builder import build_config
from neuralmonkey.attention.feed_forward import Attention
from neuralmonkey.encoders.recurrent import RecurrentEncoder
from neuralmonkey.decoders.decoder import Decoder
from neuralmonkey.decoders.encoder_projection import nematus_projection
from neuralmonkey.decoders.output_projection import nematus_output
from neuralmonkey.model.sequence import EmbeddedSequence
from neuralmonkey.vocabulary import from_nematus_json
from neuralmonkey.logging import log as _log

def log(message: str, color: str = "blue") -> None:
    _log(message, color)

def check_shape(var1_tf: tf.Variable, var2_np: np.ndarray):
    if var1_tf.get_shape().as_list() != list(var2_np.shape):
        log("Shapes do not match! Exception will follow.", color="red")

# Here come a few functions that fiddle with the Nematus parameters in order to
# fit them to Neural Monkey parameter shapes.
def emb_fix_dim1(variables: List[np.ndarray]) -> np.ndarray:
    return emb_fix(variables, dim=1)

def emb_fix(variables: List[np.ndarray], dim: int = 0) -> np.ndarray:
    """Process nematus tensors with vocabulary dimension.

    Nematus uses only two special symbols, eos and UNK. For embeddings of start
    and pad tokens, we use zero vectors, inserted to the correct position in
    the parameter matrix.

    Arguments:
        variables: the list of variables. Must be of length 1.
        dim: The vocabulary dimension.
    """
    if len(variables) != 1:
        raise ValueError("VocabFix only works with single vars. {} given."
                         .format(len(variables)))

    if dim != 0 and dim != 1:
        raise ValueError("dim can only be 0 or 1. is: {}".format(dim))

    variable = variables[0]
    shape = variable.shape

    # insert start token (hack from nematus - last from vocab - does it work? NO)
    # to_insert = np.squeeze(variable[-1] if dim == 0 else variable[:, -1])
    to_insert = np.zeros(shape[1 - dim]) if len(shape) > 1 else 0.
    variable = np.insert(variable, 0, to_insert, axis=dim)

    # insert padding token
    to_insert = np.zeros(shape[1 - dim]) if len(shape) > 1 else 0.
    variable = np.insert(variable, 0, to_insert, axis=dim)

    return variable

def sum_vars(variables: List[np.ndarray]) -> np.ndarray:
    return sum(variables)

def concat_vars(variables: List[np.ndarray]) -> np.ndarray:
    return np.concatenate(variables)

def squeeze(variables: List[np.ndarray]) -> np.ndarray:
    if len(variables) != 1:
        raise ValueError("Squeeze only works with single vars. {} given."
                         .format(len(variables)))
    return np.squeeze(variables[0])

# pylint: disable=line-too-long
# No point in line wrapping
VARIABLE_MAP = {
    "encoder_input/embedding_matrix_0": (["Wemb"], emb_fix),
    "decoder/word_embeddings": (["Wemb_dec"], emb_fix),
    "decoder/state_to_word_W": (["ff_logit_W"], emb_fix_dim1),
    "decoder/state_to_word_b": (["ff_logit_b"], emb_fix),
    "encoder/bidirectional_rnn/fw/OrthoGRUCell/gates/state_proj/kernel": (["encoder_U"], None),
    "encoder/bidirectional_rnn/fw/OrthoGRUCell/gates/input_proj/kernel": (["encoder_W"], None),
    "encoder/bidirectional_rnn/fw/OrthoGRUCell/gates/input_proj/bias": (["encoder_b"], None),
    "encoder/bidirectional_rnn/fw/OrthoGRUCell/candidate/state_proj/kernel": (["encoder_Ux"], None),
    "encoder/bidirectional_rnn/fw/OrthoGRUCell/candidate/input_proj/kernel": (["encoder_Wx"], None),
    "encoder/bidirectional_rnn/fw/OrthoGRUCell/candidate/input_proj/bias": (["encoder_bx"], None),
    "encoder/bidirectional_rnn/bw/OrthoGRUCell/gates/state_proj/kernel": (["encoder_r_U"], None),
    "encoder/bidirectional_rnn/bw/OrthoGRUCell/gates/input_proj/kernel": (["encoder_r_W"], None),
    "encoder/bidirectional_rnn/bw/OrthoGRUCell/gates/input_proj/bias": (["encoder_r_b"], None),
    "encoder/bidirectional_rnn/bw/OrthoGRUCell/candidate/state_proj/kernel": (["encoder_r_Ux"], None),
    "encoder/bidirectional_rnn/bw/OrthoGRUCell/candidate/input_proj/kernel": (["encoder_r_Wx"], None),
    "encoder/bidirectional_rnn/bw/OrthoGRUCell/candidate/input_proj/bias": (["encoder_r_bx"], None),
    "decoder/initial_state/encoders_projection/kernel": (["ff_state_W"], None),
    "decoder/initial_state/encoders_projection/bias": (["ff_state_b"], None),
    "decoder/attention_decoder/OrthoGRUCell/gates/state_proj/kernel": (["decoder_U"], None),
    "decoder/attention_decoder/OrthoGRUCell/gates/input_proj/kernel": (["decoder_W"], None),
    "decoder/attention_decoder/OrthoGRUCell/gates/input_proj/bias": (["decoder_b"], None),
    "decoder/attention_decoder/OrthoGRUCell/candidate/state_proj/kernel": (["decoder_Ux"], None),
    "decoder/attention_decoder/OrthoGRUCell/candidate/input_proj/kernel": (["decoder_Wx"], None),
    "decoder/attention_decoder/OrthoGRUCell/candidate/input_proj/bias": (["decoder_bx"], None),
    "attention/attn_key_projection": (["decoder_Wc_att"], None),
    "attention/attn_projection_bias": (["decoder_b_att"], None),
    "attention/Attention/attn_query_projection": (["decoder_W_comb_att"], None),
    "attention/attn_similarity_v": (["decoder_U_att"], squeeze),
    "attention/attn_bias": (["decoder_c_tt"], squeeze),
    "decoder/attention_decoder/cond_gru_2_cell/gates/state_proj/kernel": (["decoder_U_nl"], None),
    "decoder/attention_decoder/cond_gru_2_cell/gates/input_proj/kernel": (["decoder_Wc"], None),
    "decoder/attention_decoder/cond_gru_2_cell/gates/state_proj/bias": (["decoder_b_nl"], None),
    "decoder/attention_decoder/cond_gru_2_cell/candidate/state_proj/kernel": (["decoder_Ux_nl"], None),
    "decoder/attention_decoder/cond_gru_2_cell/candidate/input_proj/kernel": (["decoder_Wcx"], None),
    "decoder/attention_decoder/cond_gru_2_cell/candidate/state_proj/bias": (["decoder_bx_nl"], None),
    "decoder/attention_decoder/rnn_state/kernel": (["ff_logit_lstm_W"], None),
    "decoder/attention_decoder/rnn_state/bias": (["ff_logit_lstm_b"], None),
    "decoder/attention_decoder/prev_out/kernel": (["ff_logit_prev_W"], None),
    "decoder/attention_decoder/prev_out/bias": (["ff_logit_prev_b"], None),
    "decoder/attention_decoder/context/kernel": (["ff_logit_ctx_W"], None),
    "decoder/attention_decoder/context/bias": (["ff_logit_ctx_b"], None)
}
# pylint: enable=line-too-long

ENCODER_NAME = "encoder"
DECODER_NAME = "decoder"
ATTENTION_NAME = "attention"

def load_nematus_json(path: str) -> Dict:
    with open(path, "r", encoding="utf-8") as f_json:
        contents = json.load(f_json)

    prefix = os.path.realpath(os.path.dirname(path))
    config = {
        "encoder_type": contents["encoder"],
        "decoder_type": contents["decoder"],
        "n_words_src": contents["n_words_src"],
        "n_words_tgt": contents["n_words"],
        "variables_file": contents["saveto"],
        "rnn_size": contents["dim"],
        "embedding_size": contents["dim_word"],
        "src_vocabulary": os.path.join(
            prefix, contents["dictionaries"][0]),
        "tgt_vocabulary": os.path.join(
            prefix, contents["dictionaries"][1]),
        "max_length": contents["maxlen"]
    }

    if config["encoder_type"] != "gru":
        raise ValueError("Unsupported encoder type: {}"
                         .format(config["encoder_type"]))

    if config["decoder_type"] != "gru_cond":
        raise ValueError("Unsupported decoder type: {}"
                         .format(config["decoder_type"]))

    if not os.path.isfile(config["src_vocabulary"]):
        raise FileNotFoundError("Vocabulary file not found: {}"
                                .format(config["src_vocabulary"]))

    if not os.path.isfile(config["tgt_vocabulary"]):
        raise FileNotFoundError("Vocabulary file not found: {}"
                                .format(config["tgt_vocabulary"]))

    return config

VOCABULARY_TEMPLATE = """\
[vocabulary_{}]
class=vocabulary.from_nematus_json
path="{}"
max_size={}
pad_to_max_size=True
"""

ENCODER_TEMPLATE = """\
[encoder]
class=encoders.RecurrentEncoder
name="{}"
input_sequence=<input_sequence>
rnn_size={}
rnn_cell="NematusGRU"
dropout_keep_prob=1.0

[input_sequence]
class=model.sequence.EmbeddedSequence
name="{}"
vocabulary=<vocabulary_src>
data_id="source"
embedding_size={}
max_length={}
add_end_symbol=True
"""

def build_encoder(config: Dict) -> Tuple[RecurrentEncoder, str]:
    vocabulary = from_nematus_json(
        config["src_vocabulary"], max_size=config["n_words_src"],
        pad_to_max_size=True)

    vocabulary_ini = VOCABULARY_TEMPLATE.format(
        "src", config["src_vocabulary"], config["n_words_src"])

    inp_seq_name = "{}_input".format(ENCODER_NAME)
    inp_seq = EmbeddedSequence(
        name=inp_seq_name,
        vocabulary=vocabulary,
        data_id="source",
        embedding_size=config["embedding_size"])

    encoder = RecurrentEncoder(
        name=ENCODER_NAME,
        input_sequence=inp_seq,
        rnn_size=config["rnn_size"],
        rnn_cell="NematusGRU")

    encoder_ini = ENCODER_TEMPLATE.format(
        ENCODER_NAME, config["rnn_size"],
        inp_seq_name, config["embedding_size"], config["max_length"])

    return encoder, "\n".join([vocabulary_ini, encoder_ini])

ATTENTION_TEMPLATE = """\
[attention]
class=attention.Attention
name="{}"
encoder=<encoder>
dropout_keep_prob=1.0
"""

def build_attention(config: Dict,
                    encoder: RecurrentEncoder) -> Tuple[Attention, str]:
    attention = Attention(
        name=ATTENTION_NAME,
        encoder=encoder)

    attention_ini = ATTENTION_TEMPLATE.format(ATTENTION_NAME)

    return attention, attention_ini

DECODER_TEMPLATE = """\
[decoder]
class=decoders.Decoder
name="{}"
vocabulary=<vocabulary_tgt>
data_id="target"
embedding_size={}
rnn_size={}
max_output_len={}
encoders=[<encoder>]
encoder_projection=<nematus_mean>
attentions=[<attention>]
attention_on_input=False
conditional_gru=True
output_projection=<nematus_nonlinear>
rnn_cell="NematusGRU"
dropout_keep_prob=1.0

[nematus_nonlinear]
class=decoders.output_projection.nematus_output
output_size={}
dropout_keep_prob=1.0

[nematus_mean]
class=decoders.encoder_projection.nematus_projection
dropout_keep_prob=1.0
"""

def build_decoder(config: Dict,
                  attention: Attention,
                  encoder: RecurrentEncoder) -> Tuple[Decoder, str]:
    vocabulary = from_nematus_json(
        config["tgt_vocabulary"],
        max_size=config["n_words_tgt"],
        pad_to_max_size=True)

    vocabulary_ini = VOCABULARY_TEMPLATE.format(
        "tgt", config["tgt_vocabulary"], config["n_words_tgt"])

    decoder = Decoder(
        name=DECODER_NAME,
        vocabulary=vocabulary,
        data_id="target",
        max_output_len=config["max_length"],
        embedding_size=config["embedding_size"],
        rnn_size=config["rnn_size"],
        encoders=[encoder],
        attentions=[attention],
        attention_on_input=False,
        conditional_gru=True,
        encoder_projection=nematus_projection(dropout_keep_prob=1.0),
        output_projection=nematus_output(config["embedding_size"]),
        rnn_cell="NematusGRU")

    decoder_ini = DECODER_TEMPLATE.format(
        DECODER_NAME, config["embedding_size"], config["rnn_size"],
        config["max_length"], config["embedding_size"])

    return decoder, "\n".join([vocabulary_ini, decoder_ini])

def build_model(config: Dict) -> Tuple[
        RecurrentEncoder, Attention, Decoder, str]:
    encoder, encoder_cfg = build_encoder(config)
    attention, attention_cfg = build_attention(config, encoder)
    decoder, decoder_cfg = build_decoder(config, attention, encoder)

    ini = "\n".join([encoder_cfg, attention_cfg, decoder_cfg])

    return ini

def load_nematus_file(path: str) -> Dict[str, np.ndarray]:
    contents = np.load(path)
    cnt_dict = dict(contents)
    contents.close()
    return cnt_dict

def assign_vars(variables: Dict[str, np.ndarray]) -> List[tf.Tensor]:
    """For each variable in the map, assign the value from the dict"""

    trainable_vars = tf.trainable_variables()
    assign_ops = []

    for var in trainable_vars:
        map_key = var.op.name

        if map_key not in VARIABLE_MAP:
            raise ValueError("Map key {} not in variable map".format(map_key))

        nem_var_list, fun = VARIABLE_MAP[map_key]

        for nem_var in nem_var_list:
            if nem_var not in variables:
                raise ValueError("Alleged nematus var {} not found in loaded "
                                 "nematus vars.".format(nem_var))

        if fun is None:
            if len(nem_var_list) != 1:
                raise ValueError(
                    "Var list for map key {} must have length 1. "
                    "Length {} found instead."
                    .format(map_key, len(nem_var_list)))

            to_assign = variables[nem_var_list[0]]
        else:
            to_assign = fun([variables[v] for v in nem_var_list])

        check_shape(var, to_assign)
        assign_ops.append(tf.assign(var, to_assign))

    return assign_ops

INI_HEADER = """\
; This is an automatically generated configuration file
; for running imported nematus model
; For further training, set the configuration as appropriate

[main]
name="nematus imported translation"
tf_manager=<tf_manager>
output="{}"
runners=[<runner>]
postprocess=None
evaluation=[("target", evaluators.bleu.BLEU)]
runners_batch_size=1

; TODO Set these additional attributes for further training
; batch_size=80
; epochs=10
; train_dataset=<train_data>
; val_dataset=<val_data>
; trainer=<trainer>
; logging_period=20
; validation_period=60
; random_seed=1234

; [train_data]
; class=dataset.load_dataset_from_files
; s_source="PATH/TO/DATA" ; TODO do not forget to fill this out!
; s_target="PATH/TO/DATA" ; TODO do not forget to fill this out!
; lazy=True

; [val_data]
; class=dataset.load_dataset_from_files
; s_source="PATH/TO/DATA" ; TODO do not forget to fill this out!
; s_target="PATH/TO/DATA" ; TODO do not forget to fill this out!

; [trainer]
; class=trainers.cross_entropy_trainer.CrossEntropyTrainer
; decoders=[<decoder>]
; l2_weight=1.0e-8
; clip_norm=1.0

[tf_manager]
class=tf_manager.TensorFlowManager
num_threads=4
num_sessions=1

[runner]
class=runners.runner.GreedyRunner
decoder=<decoder>
output_series="target"
"""

def write_config(experiment_dir: str, ini: str) -> None:
    experiment_file = os.path.join(experiment_dir, "experiment.ini")
    with open(experiment_file, "w", encoding="utf-8") as f_out:
        f_out.write(INI_HEADER.format(experiment_dir))
        f_out.write(ini)

def prepare_output_dir(output_dir: str) -> bool:
    if os.path.isdir(output_dir):
        log("Directory {} already exists. Choose a nonexistent one.".
            format(output_dir))
        exit(1)

    os.mkdir(output_dir)

def main() -> None:
    log("Script started.")
    parser = argparse.ArgumentParser(description=__doc__)
    parser.add_argument("nematus_json", metavar="NEMATUS-JSON",
                        help="nematus json file")
    parser.add_argument("nematus_variables", metavar="NEMATUS-FILE",
                        help="nematus variable file")
    parser.add_argument("output_dir", metavar="OUTPUT-DIR",
                        help="output directory")
    args = parser.parse_args()

    log("Loading nematus variables from {}.".format(args.nematus_variables))
    nematus_vars = load_nematus_file(args.nematus_variables)

    log("Loading nematus JSON config from {}.".format(args.nematus_json))
    nematus_json_cfg = load_nematus_json(args.nematus_json)

    log("Bulding model.")
    ini = build_model(nematus_json_cfg)

    log("Defining assign Ops.")
    assign_ops = assign_vars(nematus_vars)

    log("Preparing output directory {}".format(args.output_dir))
    prepare_output_dir(args.output_dir)

    log("Writing configuration file to {}/experiment.ini."
        .format(args.output_dir))
    write_config(args.output_dir, ini)

    log("Creating TF session.")
    s = tf.Session()

    log("Running session to assign to Neural Monkey variables.")
    s.run(assign_ops)

    log("Initializing saver.")
    saver = tf.train.Saver()

    variables_file = os.path.join(args.output_dir, "variables.data")
    log("Saving variables to {}".format(variables_file))
    saver.save(s, variables_file)

    log("Finished.")

if __name__ == "__main__":
    main()

max_stars_repo_name: orangeGoran/neuralmonkey-experiments

max_issues_repo_name: orangeGoran/neuralmonkey-experiments

max_forks_repo_name: orangeGoran/neuralmonkey-experiments


text:
import os
import numpy as np
import pandas as pd
import cv2
import matplotlib.pyplot as plt

def create_mask_for_image(image):
    '''
    Utility Function to create segmented morphological masks
    '''
    image_hsv = cv2.cvtColor(image, cv2.COLOR_BGR2HSV)

    sensitivity = 35
    lower_hsv = np.array([60 - sensitivity, 100, 50])
    upper_hsv = np.array([60 + sensitivity, 255, 255])

    mask = cv2.inRange(image_hsv, lower_hsv, upper_hsv)
    kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (11,11))
    mask = cv2.morphologyEx(mask, cv2.MORPH_CLOSE, kernel)

    return mask

def segment_image(image):
    '''
    Utility Function to apply segmented morphological masks
    '''
    mask = create_mask_for_image(image)
    output = cv2.bitwise_and(image, image, mask = mask)
    return output

def sharpen_image(image):
    '''
    Utility Function to sharpen processed images
    '''
    image_blurred = cv2.GaussianBlur(image, (0, 0), 3)
    image_sharp = cv2.addWeighted(image, 1.5, image_blurred, -0.5, 0)
    return image_sharp

x_train = []

for i in range(len(train)):
    #print(train_df['file'][i])
    img = cv2.imread(train_df['file'][i])
    # cv2_imshow(img)
    # img = cv2.resize(img,dsize=(256,256))
    if i == 5 :
      print('Original Image')
      cv2_imshow(img)
    img_stack = segment_image(img)
    if i == 5 :
      print('Segmented Image')
      cv2_imshow(img_stack)
    img_stack = sharpen_image(img_stack)
    if i == 5 :
      print('Sharpened Image')
      cv2_imshow(img_stack)
    img_stack = cv2.cvtColor( img_stack, cv2.COLOR_RGB2GRAY )
    if i == 5 :
      print('Gray Image')
      cv2_imshow(img_stack)
    img_stack = np.reshape(img_stack,(1024,1024,1))
    x_train.append(np.concatenate((np.array(img),np.array(img_stack)),axis=2))

x_train = np.array(x_train)

max_stars_repo_path: Code Snippets/preprocessing.py

max_issues_repo_path: Code Snippets/preprocessing.py

max_forks_repo_path: Code Snippets/preprocessing.py


text:
import numpy as np

from nnaps.mesa.evolution_phases import _check_history_parameters, HeIgF

def mass_loss_error(history):
    """
    Check if there is a possible issue with mass loss. This function checks if RLOF is still ongoing when the model
    stops. For a stable model this would be an indication that something is wrong.

    Required history parameters:
        - lg_mstar_dot_1
        - lg_wind_mdot_1

    :param history: the evolution history of the system.
    :type history: numpy ndarray
    :return: True if there is a mass loss error, else False
    """

    required_parameters = ['lg_mstar_dot_1', 'lg_wind_mdot_1']
    if not _check_history_parameters(history, required_parameters, evol_phase='ML error', raise_error=False):
        return False

    mass_loss = np.log10(10 ** history['lg_mstar_dot_1'] - 10 ** history['lg_wind_mdot_1'])

    return mass_loss[-1] > -10

def he_ignition_error(history):
    """
    Check if there are issues with He ignition. The function checks if there is ignition but no actual core He
    burning, and also if there is a ramp up to ignition, but the ignition itself is unsuccessful.

    Required history parameters:
        - age
        - log_LHe
        - log_center_T
        - log_center_Rho

    :param history: the evolution history of the system.
    :type history: numpy ndarray
    :return: True if there is a He ignition error, else False
    """
    required_parameters = ['age', 'log_LHe', 'log_center_T', 'log_center_Rho']
    if not _check_history_parameters(history, required_parameters, evol_phase='He ignition error', raise_error=False):
        return False

    # Check if there is a ramp-up to He ignition, but no actual ignition. Check the last 2000 points
    if np.all(history['log_center_T'] < HeIgF(history['log_center_Rho'])) and np.all(history['log_LHe'] < 1):
        d = history[history['age'] > history[history['log_LHe'] > -50]['age'][0]]
        if len(d) > 5000:
            log_LHe = d['log_LHe'][-5000:]
            log_LHe_avg = np.average(log_LHe)
            if (log_LHe > log_LHe_avg - 0.5).all() and (log_LHe < log_LHe_avg + 0.5).all() \
                    and log_LHe_avg > -50 and log_LHe_avg > history['log_LHe'][0]:
                return True

    # Check if there is He ignition but no core burning
    if np.any(history['log_LHe'] > 1) and np.all(history['log_center_T'] < HeIgF(history['log_center_Rho'])):
        return True

    return False

def he_core_burning_error(history):
    """
    Check if there is an issue with He core burning. The function checks if core He burning starts, but never ends
    before the end of the evolution.

    Required history parameters:
        - c_core_mass
        - log_center_T
        - log_center_Rho

    :param history: the evolution history of the system.
    :type history: numpy ndarray
    :return: True if there is a He core burning error, else False
    """
    required_parameters = ['c_core_mass', 'log_center_T', 'log_center_Rho']
    if not _check_history_parameters(history, required_parameters, evol_phase='He ignition error', raise_error=False):
        return False

    if np.any(history['log_center_T'] >= HeIgF(history['log_center_Rho'])) and np.all(history['c_core_mass'] < 0.01):
        return True
    else:
        return False

def check_error_flags(history, termination_code):
    """
    Check for some possible errors in the model and report them.

    Errors that are checked and there codes:

        1. MESA stopped because of max model number
        2. MESA stopped because the accretor is overflowing it's Roche-lobe
        3. The mass loss phase doesn't end, to check if mass loss didn't cause mesa to crash:
           :func:`~nnaps.mesa.evolution_errors.mass_loss_error`
        4. Potential problem with He ignition: He ramps up and doesnt ignite, or He ignites but there is no core
           burning: :func:`~nnaps.mesa.evolution_errors.he_ignition_error`
        5. He core burning starts, but there is no formation of a CO core:
           :func:`~nnaps.mesa.evolution_errors.he_core_burning_error`

    :param history: the evolution history of the system
    :type history: numpy ndarray
    :param termination_code: the termination code of the mesa model
    :type termination_code: str
    :return: list of integers with error codes
    """

    error_codes = []

    # Check termination_code
    if termination_code == 'max_model_number':
        error_codes.append(1)
    if termination_code == 'accretor_overflow_terminate':
        error_codes.append(2)

    # Check Mass loss error
    if mass_loss_error(history):
        error_codes.append(3)

    # Check He ignition error
    if he_ignition_error(history):
        error_codes.append(4)

    # check He core burning errors
    if he_core_burning_error(history):
        error_codes.append(5)

    return error_codes

max_stars_repo_name: vosjo/nnaps

max_issues_repo_name: vosjo/nnaps

max_forks_repo_name: vosjo/nnaps