text:
from astropy.coordinates import Angle, SkyCoord
from regions import CircleSkyRegion
import matplotlib.pyplot as plt
from gammapy.makers import ReflectedRegionsFinder
from gammapy.maps import RegionGeom, WcsNDMap

# Exclude a rectangular region
exclusion_mask = WcsNDMap.create(npix=(801, 701), binsz=0.01, skydir=(83.6, 23.0))

coords = exclusion_mask.geom.get_coord().skycoord
data = (Angle("23 deg") < coords.dec) & (coords.dec < Angle("24 deg"))
exclusion_mask.data = ~data

pos = SkyCoord(83.633, 22.014, unit="deg")
radius = Angle(0.3, "deg")
on_region = CircleSkyRegion(pos, radius)
center = SkyCoord(83.633, 24, unit="deg")

# One can impose a minimal distance between ON region and first reflected regions
finder = ReflectedRegionsFinder(
    min_distance_input="0.2 rad",
)
regions, _ = finder.run(
    region=on_region,
    center=center,
    exclusion_mask=exclusion_mask,
)

fig, axes = plt.subplots(
    ncols=3,
    subplot_kw={"projection": exclusion_mask.geom.wcs},
    figsize=(12, 3),
)

def plot_regions(ax, regions, on_region, exclusion_mask):
    """Little helper function to plot off regions"""
    exclusion_mask.plot_mask(ax=ax, colors="gray")
    on_region.to_pixel(ax.wcs).plot(ax=ax, edgecolor="tab:orange")
    geom = RegionGeom.from_regions(regions)
    geom.plot_region(ax=ax, color="tab:blue")

ax = axes[0]
ax.set_title("Min. distance first region")
plot_regions(ax=ax, regions=regions, on_region=on_region, exclusion_mask=exclusion_mask)

# One can impose a minimal distance between two adjacent regions
finder = ReflectedRegionsFinder(
    min_distance="0.1 rad",
)
regions, _ = finder.run(
    region=on_region,
    center=center,
    exclusion_mask=exclusion_mask,
)

ax = axes[1]
ax.set_title("Min. distance all regions")
plot_regions(ax=ax, regions=regions, on_region=on_region, exclusion_mask=exclusion_mask)

# One can impose a maximal number of regions to be extracted
finder = ReflectedRegionsFinder(
    max_region_number=5,
    min_distance="0.1 rad",
)
regions, _ = finder.run(
    region=on_region,
    center=center,
    exclusion_mask=exclusion_mask,
)

ax = axes[2]
ax.set_title("Max. number of regions")
plot_regions(ax=ax, regions=regions, on_region=on_region, exclusion_mask=exclusion_mask)
plt.show()

max_stars_repo_name: aaguasca/gammapy

max_issues_repo_name: aaguasca/gammapy

max_forks_repo_name: aaguasca/gammapy

reason: from astropy


text:
import matplotlib.pyplot as plt
import numpy as np
from environment.gride_world import GridWorld

import TebularPlanningAndLearning.Dyna_Q_plus as DQp

if __name__ == '__main__':
    episode_numbers = 400
    dqp_reward = 0.0
    average_times = 10
    average_step_array = np.zeros([2, episode_numbers])
    for ave_iter in range(average_times):
        env = GridWorld(6, [25, 26, 27, 28, 29], start_position=32, end_position_list=[5])
        agent_dqp = DQp.Agent(env, n=10, kappa=0.001)
        agent_dq = DQp.Agent(env, n=10, kappa=0)
        dqp_step_rewards_list = []
        dq_step_rewards_list = []
        for i in range(episode_numbers):
            agent_dqp.dyna_q_plus(1, alpha=0.1, gamma=0.95, epsilon=.3)
            dqp_step_rewards_list.append(agent_dqp.total_step_num)
            agent_dq.dyna_q_plus(1, alpha=0.1, gamma=0.95, epsilon=.3)
            dq_step_rewards_list.append(agent_dq.total_step_num)
            if i == episode_numbers / 2:
                agent_dqp.env = GridWorld(6, [25, 26, 27, 28], start_position=32, end_position_list=[5])
                agent_dq.env = GridWorld(6, [25, 26, 27, 28], start_position=32, end_position_list=[5])
        average_step_array[0] += np.array(dqp_step_rewards_list)
        average_step_array[1] += np.array(dq_step_rewards_list)
    plt.plot((average_step_array[0] - average_step_array[0][3]) / average_times,
             np.arange(1, episode_numbers + 1), label='Dyna_Q+')

    plt.plot((average_step_array[1] - average_step_array[1][3]) / average_times,
             np.arange(1, episode_numbers + 1), label='Dyna_Q')
    plt.legend()
    plt.show()

ext: py

reason: import numpy


text:
"""
Plotly backend for plotting methods
===================================

"""
import functools
from collections.abc import Iterable

import numpy
import pandas as pd
import plotly.express as px
import plotly.graph_objects as go
from plotly.subplots import make_subplots

import orion.analysis
import orion.analysis.base
from orion.algo.space import Categorical, Fidelity
from orion.core.worker.transformer import build_required_space

def lpi(
    experiment,
    with_evc_tree=True,
    model="RandomForestRegressor",
    model_kwargs=None,
    n_points=20,
    n_runs=10,
    **kwargs,
):
    """Plotly implementation of `orion.plotting.lpi`"""
    if not experiment:
        raise ValueError("Parameter 'experiment' is None")

    if model_kwargs is None:
        model_kwargs = {}

    df = experiment.to_pandas(with_evc_tree=with_evc_tree)
    df = df.loc[df["status"] == "completed"]

    if df.empty:
        return go.Figure()

    df = orion.analysis.lpi(
        df,
        experiment.space,
        model=model,
        n_points=n_points,
        n_runs=n_runs,
        **model_kwargs,
    )

    fig = go.Figure(
        data=[
            go.Bar(
                x=df.index.tolist(),
                y=df["LPI"].tolist(),
                error_y=dict(type="data", array=df["STD"].tolist()),
            )
        ]
    )

    y_axis_label = "Local Parameter Importance (LPI)"
    fig.update_layout(
        title=f"LPI for experiment '{experiment.name}'",
        xaxis_title="Hyperparameters",
        yaxis_title=y_axis_label,
    )

    return fig

def parallel_coordinates(
    experiment, with_evc_tree=True, order=None, colorscale="YlOrRd", **kwargs
):
    """Plotly implementation of `orion.plotting.parallel_coordinates`"""

    def build_frame():
        """Builds the dataframe for the plot"""
        names = list(experiment.space.keys())

        df = experiment.to_pandas(with_evc_tree=with_evc_tree)
        df = df.loc[df["status"] == "completed"]

        if df.empty:
            return df

        df[names] = df[names].transform(
            functools.partial(_curate_params, space=experiment.space)
        )

        df = _flatten_dims(df, experiment.space)

        return df

    def infer_order(space, order):
        """Create order if not passed, otherwise verify it"""
        params = orion.analysis.base.flatten_params(space, order)
        if order is None:
            fidelity_dims = [
                dim for dim in experiment.space.values() if isinstance(dim, Fidelity)
            ]
            fidelity = fidelity_dims[0].name if fidelity_dims else None
            if fidelity in params:
                del params[params.index(fidelity)]
                params.insert(0, fidelity)

        return params

    def get_dimension(data, name, dim):
        dim_data = dict(label=name, values=data[name])
        if dim.type == "categorical":
            categories = dim.interval()
            dim_data["tickvals"] = list(range(len(categories)))
            dim_data["ticktext"] = categories
        else:
            dim_data["range"] = dim.interval()
        return dim_data

    if not experiment:
        raise ValueError("Parameter 'experiment' is None")

    df = build_frame()

    if df.empty:
        return go.Figure()

    trial = experiment.fetch_trials_by_status("completed")[0]

    flattened_space = build_required_space(
        experiment.space, shape_requirement="flattened"
    )

    dimensions = [
        get_dimension(df, name, flattened_space[name])
        for name in infer_order(experiment.space, order)
    ]

    objective_name = trial.objective.name

    objectives = df["objective"]
    omin = min(df["objective"])
    omax = max(df["objective"])

    dimensions.append(dict(label=objective_name, range=(omin, omax), values=objectives))

    fig = go.Figure(
        data=go.Parcoords(
            line=dict(
                color=objectives,
                colorscale=colorscale,
                showscale=True,
                cmin=omin,
                cmax=omax,
                colorbar=dict(title=objective_name),
            ),
            dimensions=dimensions,
        )
    )

    fig.update_layout(
        title=f"Parallel Coordinates Plot for experiment '{experiment.name}'"
    )

    return fig

def rankings(experiments, with_evc_tree=True, order_by="suggested", **kwargs):
    """Plotly implementation of `orion.plotting.rankings`"""

    def reformat_competitions(experiments):
        if isinstance(experiments, dict) and isinstance(
            next(iter(experiments.values())), Iterable
        ):
            competitions = []
            remaining = True
            n_competitions = len(next(iter(experiments.values())))
            for ith_competition in range(n_competitions):
                competition = {}
                for name in experiments.keys():
                    competition[name] = experiments[name][ith_competition]
                competitions.append(competition)
        elif isinstance(experiments, dict):
            competitions = experiments
        elif isinstance(experiments, Iterable) and not isinstance(experiments[0], dict):
            competitions = {
                f"{experiment.name}-v{experiment.version}": experiment
                for experiment in experiments
            }
        else:
            competitions = experiments

        return competitions

    def build_groups(competitions):

        if not isinstance(competitions, dict):
            rankings = []
            for competition in competitions:
                rankings.append(build_frame(competition))
            df = pd.concat(rankings)
            data_frames = orion.analysis.average(
                df, group_by=["order", "name"], key="rank", return_var=True
            )
        else:
            data_frames = build_frame(competitions)

        return data_frames

    def build_frame(competition):
        """Builds the dataframe for the plot"""

        frames = []
        for name, experiment in competition.items():
            df = experiment.to_pandas(with_evc_tree=with_evc_tree)
            df = df.loc[df["status"] == "completed"]
            df = df.sort_values(order_by)
            df = orion.analysis.regret(df)
            df["name"] = [name] * len(df)
            df["order"] = range(len(df))
            frames.append(df)

        df = pd.concat(frames)

        return orion.analysis.ranking(df)

    def get_objective_name(competition):
        """Infer name of objective based on trials of one experiment"""
        if not isinstance(competition, dict):
            return get_objective_name(competition[0])

        for experiment in competition.values():
            trials = experiment.fetch_trials_by_status("completed")
            if trials:
                return trials[0].objective.name
        return "objective"

    ORDER_KEYS = ["suggested", "reserved", "completed"]

    if not experiments:
        raise ValueError("Parameter 'experiment' is None")

    if order_by not in ORDER_KEYS:
        raise ValueError(f"Parameter 'order_by' is not one of {ORDER_KEYS}")

    competitions = reformat_competitions(experiments)
    df = build_groups(competitions)

    fig = go.Figure()

    if df.empty:
        return fig

    names = set(df["name"])
    for i, name in enumerate(sorted(names)):
        exp_data = df[df["name"] == name]
        if "rank_mean" in exp_data:
            y = exp_data["rank_mean"]
        else:
            y = exp_data["rank"]
        x = list(range(len(y)))
        fig.add_scatter(
            x=x,
            y=y,
            mode="lines",
            line=dict(color=px.colors.qualitative.G10[i]),
            name=name,
        )
        if "rank_var" in exp_data:
            dy = exp_data["rank_var"]
            fig.add_scatter(
                x=list(x) + list(x)[::-1],
                y=list(y - dy) + list(y + dy)[::-1],
                fill="toself",
                showlegend=False,
                line=dict(color=px.colors.qualitative.G10[i], width=0),
                name=name,
            )

    objective = get_objective_name(competitions)
    fig.update_layout(
        title=f"Average Rankings",
        xaxis_title=f"Trials ordered by {order_by} time",
        yaxis_title=f"Ranking based on {objective}",
        hovermode="x",
    )
    return fig

def partial_dependencies(
    experiment,
    with_evc_tree=True,
    params=None,
    smoothing=0.85,
    n_grid_points=10,
    n_samples=50,
    colorscale="Blues",
    model="RandomForestRegressor",
    model_kwargs=None,
    verbose_hover=True,
):
    """Plotly implementation of `orion.plotting.partial_dependencies`"""

    def build_data():
        """Builds the dataframe for the plot"""
        df = experiment.to_pandas(with_evc_tree=with_evc_tree)

        names = list(experiment.space.keys())
        df["params"] = df[names].apply(_format_hyperparameters, args=(names,), axis=1)

        df = df.loc[df["status"] == "completed"]
        data = orion.analysis.partial_dependency(
            df,
            experiment.space,
            params=params,
            model=model,
            n_grid_points=n_grid_points,
            n_samples=n_samples,
            **model_kwargs,
        )
        df = _flatten_dims(df, experiment.space)
        return (df, data)

    def _set_scale(figure, dims, x, y):
        for axis, dim in zip("xy", dims):
            if "reciprocal" in dim.prior_name or dim.type == "fidelity":
                getattr(figure, f"update_{axis}axes")(type="log", row=y, col=x)

    def _plot_marginalized_avg(data, x_name):
        return go.Scatter(
            x=data[0][x_name],
            y=data[1],
            mode="lines",
            name=None,
            showlegend=False,
            line=dict(
                color=px.colors.qualitative.D3[0],
            ),
        )

    def _plot_marginalized_std(data, x_name):
        return go.Scatter(
            x=list(data[0][x_name]) + list(data[0][x_name])[::-1],
            y=list(data[1] - data[2]) + list((data[1] + data[2]))[::-1],
            mode="lines",
            name=None,
            fill="toself",
            showlegend=False,
            line=dict(
                color=px.colors.qualitative.D3[0],
                width=0,
            ),
        )

    def _plot_contour(data, x_name, y_name):
        return go.Contour(
            x=data[0][x_name],
            y=data[0][y_name],
            z=data[1],
            connectgaps=True,
            # Share the same color range across contour plots
            coloraxis="coloraxis",
            line_smoothing=smoothing,
            # To show labels
            contours=dict(
                coloring="heatmap",
                showlabels=True,  # show labels on contours
                labelfont=dict(  # label font properties
                    size=12,
                    color="white",
                ),
            ),
        )

    def _plot_scatter(x, y, df):
        return go.Scatter(
            x=x,
            y=y,
            marker={
                "line": {"width": 0.5, "color": "Grey"},
                "color": "black",
                "size": 5,
            },
            mode="markers",
            opacity=0.5,
            showlegend=False,
            customdata=list(zip(df["id"], df["suggested"], df["params"])),
            hovertemplate=_template_trials(verbose_hover),
        )

    if model_kwargs is None:
        model_kwargs = {}

    df, data = build_data()

    if not data:
        return go.Figure()

    params = [
        param_names for param_names in data.keys() if isinstance(param_names, str)
    ]

    flattened_space = build_required_space(
        experiment.space,
        shape_requirement="flattened",
    )

    fig = make_subplots(
        rows=len(params),
        cols=len(params),
        shared_xaxes=True,
        shared_yaxes=False,
    )

    fig.update_layout(paper_bgcolor="rgba(0,0,0,0)", plot_bgcolor="rgba(0,0,0,0)")
    cmin = float("inf")
    cmax = -float("inf")

    for x_i in range(len(params)):
        x_name = params[x_i]
        fig.add_trace(
            _plot_marginalized_avg(data[x_name], x_name),
            row=x_i + 1,
            col=x_i + 1,
        )
        fig.add_trace(
            _plot_marginalized_std(data[x_name], x_name),
            row=x_i + 1,
            col=x_i + 1,
        )
        fig.add_trace(
            _plot_scatter(df[x_name], df["objective"], df),
            row=x_i + 1,
            col=x_i + 1,
        )

        _set_scale(fig, [flattened_space[x_name]], x_i + 1, x_i + 1)

        fig.update_xaxes(title_text=x_name, row=len(params), col=x_i + 1)
        if x_i > 0:
            fig.update_yaxes(title_text=x_name, row=x_i + 1, col=1)
        else:
            fig.update_yaxes(title_text="Objective", row=x_i + 1, col=x_i + 1)

        for y_i in range(x_i + 1, len(params)):
            y_name = params[y_i]
            fig.add_trace(
                _plot_contour(
                    data[(x_name, y_name)],
                    x_name,
                    y_name,
                ),
                row=y_i + 1,
                col=x_i + 1,
            )
            fig.add_trace(
                _plot_scatter(df[x_name], df[y_name], df),
                row=y_i + 1,
                col=x_i + 1,
            )

            cmin = min(cmin, data[(x_name, y_name)][1].min())
            cmax = max(cmax, data[(x_name, y_name)][1].max())

            _set_scale(
                fig,
                [flattened_space[name] for name in [x_name, y_name]],
                x_i + 1,
                y_i + 1,
            )

    for x_i in range(len(params)):
        plot_id = len(params) * x_i + x_i + 1
        if plot_id > 1:
            key = f"yaxis{plot_id}_range"
        else:
            key = "yaxis_range"
        fig.update_layout(**{key: [cmin, cmax]})

    fig.update_layout(
        title=f"Partial dependencies for experiment '{experiment.name}'",
    )
    fig.layout.coloraxis.colorbar.title = "Objective"

    fig.update_layout(coloraxis=dict(colorscale=colorscale), showlegend=False)

    return fig

def regret(
    experiment, with_evc_tree=True, order_by="suggested", verbose_hover=True, **kwargs
):
    """Plotly implementation of `orion.plotting.regret`"""

    def build_frame():
        """Builds the dataframe for the plot"""
        df = experiment.to_pandas(with_evc_tree=with_evc_tree)

        names = list(experiment.space.keys())
        df["params"] = df[names].apply(_format_hyperparameters, args=(names,), axis=1)

        df = df.loc[df["status"] == "completed"]
        df = df.sort_values(order_by)
        df = orion.analysis.regret(df)
        return df

    ORDER_KEYS = ["suggested", "reserved", "completed"]

    if not experiment:
        raise ValueError("Parameter 'experiment' is None")

    if order_by not in ORDER_KEYS:
        raise ValueError(f"Parameter 'order_by' is not one of {ORDER_KEYS}")

    df = build_frame()

    fig = go.Figure()

    if df.empty:
        return fig

    trial = experiment.fetch_trials_by_status("completed")[0]

    fig.add_scatter(
        y=df["objective"],
        mode="markers",
        name="trials",
        customdata=list(zip(df["id"], df[order_by], df["params"])),
        hovertemplate=_template_trials(verbose_hover),
    )
    fig.add_scatter(
        y=df["best"],
        mode="lines",
        name="best-to-date",
        customdata=list(zip(df["best_id"], df["best"])),
        hovertemplate=_template_best(),
    )

    if trial is None:
        y_axis_label = "Objective unknown"
    else:
        y_axis_label = f"{trial.objective.type.capitalize()} '{trial.objective.name}'"

    fig.update_layout(
        title=f"Regret for experiment '{experiment.name}'",
        xaxis_title=f"Trials ordered by {order_by} time",
        yaxis_title=y_axis_label,
    )

    return fig

def regrets(experiments, with_evc_tree=True, order_by="suggested", **kwargs):
    """Plotly implementation of `orion.plotting.regrets`"""

    compute_average = bool(
        isinstance(experiments, dict)
        and isinstance(next(iter(experiments.values())), Iterable)
    )

    def build_groups():
        """Build dataframes for groups of experiments"""
        # TODO move this
        if compute_average:
            data_frames = dict()
            for name, group in experiments.items():
                df = orion.analysis.average(build_frame(group), return_var=True)
                df["name"] = [name] * len(df)
                data_frames[name] = df
            data_frames = pd.concat(data_frames)
        elif isinstance(experiments, dict):
            data_frames = build_frame(experiments.values(), list(experiments.keys()))
        else:
            data_frames = build_frame(experiments)

        return data_frames

    def build_frame(experiments, names=None):
        """Builds the dataframe for the plot"""
        frames = []
        for i, experiment in enumerate(experiments):
            df = experiment.to_pandas(with_evc_tree=with_evc_tree)
            df = df.loc[df["status"] == "completed"]
            df = df.sort_values(order_by)
            df = orion.analysis.regret(df)
            df["name"] = [
                names[i] if names else f"{experiment.name}-v{experiment.version}"
            ] * len(df)
            df["order"] = range(len(df))
            frames.append(df)

        return pd.concat(frames)

    def get_objective_name(experiments):
        """Infer name of objective based on trials of one experiment"""
        if compute_average and isinstance(experiments, dict):
            return get_objective_name(sum(map(list, experiments.values()), []))

        if isinstance(experiments, dict):
            experiments = experiments.values()

        for experiment in experiments:
            trials = experiment.fetch_trials_by_status("completed")
            if trials:
                return trials[0].objective.name
        return "objective"

    ORDER_KEYS = ["suggested", "reserved", "completed"]

    if not experiments:
        raise ValueError("Parameter 'experiment' is None")

    if order_by not in ORDER_KEYS:
        raise ValueError(f"Parameter 'order_by' is not one of {ORDER_KEYS}")

    df = build_groups()

    fig = go.Figure()

    if df.empty:
        return fig

    names = set(df["name"])
    for i, name in enumerate(sorted(names)):
        exp_data = df[df["name"] == name]
        if "best_mean" in exp_data:
            y = exp_data["best_mean"]
        else:
            y = exp_data["best"]
        x = list(range(len(y)))
        fig.add_scatter(
            x=x,
            y=y,
            mode="lines",
            line=dict(color=px.colors.qualitative.G10[i]),
            name=name,
            legendgroup=name,
        )
        if "best_var" in exp_data:
            dy = numpy.sqrt(exp_data["best_var"])
            fig.add_scatter(
                x=list(x) + list(x)[::-1],
                y=list(y - dy) + list(y + dy)[::-1],
                fill="toself",
                showlegend=False,
                line=dict(color=px.colors.qualitative.G10[i]),
                name=name,
                legendgroup=name,
            )

    fig.update_layout(
        title=f"Average Regret",
        xaxis_title=f"Trials ordered by {order_by} time",
        yaxis_title=get_objective_name(experiments),
        hovermode="x",
    )

    return fig

def _format_value(value):
    """
    Hyperparameter can have many types, sometimes they can even be lists.
    If one of the value is a float, it has to be compact.
    """
    if isinstance(value, str):
        return value

    try:
        return f"[{','.join(_format_value(x) for x in value)}]"

    except TypeError:
        if isinstance(value, float):
            return f"{value:.5G}"
        else:
            return value

def _format_hyperparameters(hyperparameters, names):
    result = ""

    for name, value in zip(names, hyperparameters):
        x = f"<br>  {name}: {_format_value(value)}"
        result += x

    return result

def _template_trials(verbose_hover):
    template = (
        "<b>ID: %{customdata[0]}</b><br>"
        "value: %{y}<br>"
        "time: %{customdata[1]|%Y-%m-%d %H:%M:%S}<br>"
    )

    if verbose_hover:
        template += "parameters: %{customdata[2]}"

    template += "<extra></extra>"

    return template

def _template_best():
    return (
        "<b>Best ID: %{customdata[0]}</b><br>"
        "value: %{customdata[1]}"
        "<extra></extra>"
    )

def _curate_params(data, space):
    dim = space[data.name]
    if isinstance(dim, Categorical):
        data = numpy.array(data.tolist())  # To unpack lists if dim shape > 1
        shape = data.shape
        assert len(shape) <= 2
        idx = numpy.argmax(data.reshape(-1, 1) == numpy.array(dim.categories), axis=1)
        idx = idx.reshape(shape)
        if len(shape) > 1:
            return [list(idx[i]) for i in range(shape[0])]
        return idx
    return data

def _flatten_dims(data, space):
    for key, dim in space.items():
        if dim.shape:
            assert len(dim.shape) == 1

            # expand df.tags into its own dataframe
            values = data[key].apply(pd.Series)

            # rename each hp
            values = values.rename(columns=lambda x: f"{key}[{x}]")

            data = pd.concat([data[:], values[:]], axis=1)

    return data

max_stars_repo_name: nurbal/orion

max_issues_repo_name: nurbal/orion

max_forks_repo_name: mirkobronzi/orion


text:
from wows_shell import *
import numpy as np
import matplotlib
#matplotlib.use('TkAgg',warn=False, force=True)
import matplotlib.pyplot as plt

s = shell(shellParams(.460, 780, .292, 1460,
                      2574, 6, .033, 76, 45, 60, 0), "Yamato")
impacts = {}

c = shellCalc()
c.setDtMin(.01)
c.calcImpactForwardEuler(s)
impacts['Forward Euler .01'] = s.getImpact()

c.setDtMin(.04)
c.calcImpactRungeKutta2(s)
impacts['Runge Kutta .01'] = s.getImpact()

c.setDtMin(.0001)
c.calcImpactRungeKutta4(s)
reference = s.getImpact()

c.calcImpactForwardEuler(s)
impacts['Forward Euler .0001'] = s.getImpact()

c.calcImpactRungeKutta2(s)
impacts['RungeKutta2 .0001'] = s.getImpact()

plt.plot(impacts['Forward Euler .01'][impactIndices.distance, :],
         impacts['Forward Euler .01'][impactIndices.effectivePenetrationHorizontalNormalized, :])
plt.plot(impacts['Runge Kutta .01'][impactIndices.distance, :],
         impacts['Runge Kutta .01'][impactIndices.effectivePenetrationHorizontalNormalized, :])
plt.plot(reference[impactIndices.distance, :],
         reference[impactIndices.effectivePenetrationHorizontalNormalized, :])

plt.show()


text:
import os
from os.path import join
import re
import numpy as np
import time as time
from scipy.io import FortranFile
import ctypes as c
import struct

class ReadTreebrick_lowp:
    #low precision Treebricks (mostly HAGN)
    def __init__(self,file_path=None):

        """
        Make a treebricks dictionary out of file
        """
        self.file_path = file_path
        self.treebricks_dict = None
        self.read_data()

    def read_data(self):
        t0 = time.time()
        f = FortranFile(self.file_path, 'r')

        nbodies = f.read_record('i')
        mpart = f.read_record('f')
        aexp = f.read_record('f')
        omega_t = f.read_record('f')
        age_univ = f.read_record('f')
        nsub = f.read_record('i')

        self.treebricks_dict = {}
        self.treebricks_dict['nbodies'] = nbodies
        self.treebricks_dict['mpart'] = mpart
        self.treebricks_dict['aexp'] = aexp
        self.treebricks_dict['omega_t'] = omega_t
        self.treebricks_dict['age_univ'] = age_univ
        self.treebricks_dict['nh_old'] = nsub[0]
        self.treebricks_dict['nsub_old'] = nsub[1]
        nhaloes = nsub[0]+nsub[1]
        self.treebricks_dict['nhaloes'] = nhaloes
        #initialize empty list to hold all halos
        self.treebricks_dict['haloes'] = []

        print('nbodies:',nbodies,
              'mpart:', mpart,
              'aexp:', aexp,
              'omega_t:', omega_t,
              'age:', age_univ,
              'nsub:',nsub,
              'nhaloes:', nhaloes)

        nb_of_haloes = self.treebricks_dict['nh_old']+self.treebricks_dict['nsub_old']

        #define the length of the box in Mpc for New Horizon
        lbox_NH=2.0*0.07*(142.8571428)*aexp
        lbox_HAGN = (142.8571428)*aexp
        halfbox=lbox_NH/2.0
        self.treebricks_dict['lbox_NH'] = lbox_NH
        self.treebricks_dict['lbox_HAGN'] = lbox_HAGN

        #def read_halo():
        t1 = time.time()
        for i in range(nhaloes):
            halo_dict = {}

            npart = f.read_record('i')
            halo_dict['npart'] = npart.tolist()

            #members = f.read_record('i').reshape(npart, order='F')
            members = f.read_record('i')
            halo_dict['members'] = members.tolist()

            my_number = f.read_record('i')
            halo_dict['my_number'] = my_number.tolist()

            my_timestep = f.read_record('i')
            halo_dict['my_timestep'] = my_timestep.tolist()

            level_ids = f.read_record('i')
            level_ids = level_ids.tolist()

            level,host_halo,host_sub,nchild,nextsub = level_ids[0],level_ids[1],level_ids[2],level_ids[3],level_ids[4]
            halo_dict['level'] = level
            halo_dict['host_halo']  = host_halo
            halo_dict['host_sub'] = host_sub
            halo_dict['nchild'] = nchild
            halo_dict['nextsub'] = nextsub

            mass = f.read_record('f') #d for NH
            halo_dict['mass'] = mass.tolist()

            p = f.read_record('f')
            p = p.tolist()
            py,px,pz = p[0],p[1],p[2]
            halo_dict['px'] = px
            halo_dict['py'] = py
            halo_dict['pz'] = pz

            v = f.read_record('f')
            v = v.tolist()
            vx,vy,vz = v[0],v[1],v[2]
            halo_dict['vx'] = vx
            halo_dict['vy'] = vy
            halo_dict['vz'] = vz

            L = f.read_record('f')
            L = L.tolist()
            Lx,Ly,Lz = L[0],L[1],L[2]
            halo_dict['Lx'] = Lx
            halo_dict['Ly'] = Ly
            halo_dict['Lz'] = Lz

            shape = f.read_record('f')
            shape = shape.tolist()
            rmax,a,b,c = shape[0],shape[1],shape[2],shape[3]
            halo_dict['rmax'] = rmax
            halo_dict['a'] = a
            halo_dict['b'] = b
            halo_dict['c'] = c

            energy = f.read_record('f')
            energy = energy.tolist()
            ek,ep,et = energy[0],energy[1],energy[2]
            halo_dict['ek'] = ek
            halo_dict['ep'] = ep
            halo_dict['et'] = et

            spin = f.read_record('f')
            halo_dict['spin'] = spin.tolist()

            virial = f.read_record('f')
            virial = virial.tolist()
            rvir,mvir,tvir,cvel = virial[0],virial[1],virial[2],virial[3]
            halo_dict['rvir'] = rvir
            halo_dict['mvir'] = mvir
            halo_dict['tvir'] = tvir
            halo_dict['cvel'] = cvel

            halo_profile = f.read_record('f')
            halo_profile = halo_profile.tolist()
            rho_0, r_c = halo_profile[0],halo_profile[1]
            halo_dict['rho_0'] = rho_0
            halo_dict['r_c'] = r_c

            #Positions are in Mpc, we now put them back in "code units", that is assuming the length of the simulation is 1.
            #for New Horizon scale
            """

            px_NH = (halo_dict['px']/self.treebricks_dict['lbox_NH']) + 0.5
            halo_dict['px_NH'] = px_NH

            py_NH = (halo_dict['py']/self.treebricks_dict['lbox_NH']) + 0.5
            halo_dict['py_NH'] = py_NH

            pz_NH = (halo_dict['pz']/self.treebricks_dict['lbox_NH']) + 0.5
            halo_dict['pz_NH'] = pz_NH

            rvir_NH = (halo_dict['rvir']/self.treebricks_dict['lbox_NH'])
            halo_dict['rvir_NH'] = rvir_NH

            #for Horizon AGN scale
            px_HAGN = (halo_dict['px']/self.treebricks_dict['lbox_HAGN']) + 0.5
            halo_dict['px_HAGN'] = px_HAGN

            py_HAGN = (halo_dict['py']/self.treebricks_dict['lbox_HAGN']) + 0.5
            halo_dict['py_HAGN'] = py_HAGN

            pz_HAGN = (halo_dict['pz']/self.treebricks_dict['lbox_HAGN']) + 0.5
            halo_dict['pz_HAGN'] = pz_HAGN

            rvir_HAGN = (halo_dict['rvir']/self.treebricks_dict['lbox_HAGN'])
            halo_dict['rvir_HAGN'] = rvir_HAGN
            """

            #add halo dict to main dict under key of my_number

            #return halo_dict
            self.treebricks_dict['haloes'].append(halo_dict)

        #t1 = time.time()
        #for i in range(nhaloes):
            #halo_dict = read_halo()
            #self.treebricks_dict['haloes'].append(halo_dict)

        t2 = time.time()
        print('Reading haloes took {:0.2f} secs.'.format(t2-t1))
        print('Total time was {:0.2f} secs.'.format(t2-t0))

        return self.treebricks_dict

class GalaxyCatalog:
    #High precision Treebricks (NH)
    def __init__(self,file_path=None):

        """
        Make a treebricks dictionary out of file
        """
        self.file_path = file_path
        self.treebricks_dict = None
        self.read_data()

    def read_data(self):
        t0 = time.time()
        f = FortranFile(self.file_path, 'r')

        nbodies = f.read_record('i')
        mpart = f.read_record('d')
        aexp = f.read_record('d')
        omega_t = f.read_record('d')
        age_univ = f.read_record('d')
        nsub = f.read_record('i')

        self.treebricks_dict = {}
        self.treebricks_dict['nbodies'] = nbodies
        self.treebricks_dict['mpart'] = mpart
        self.treebricks_dict['aexp'] = aexp
        self.treebricks_dict['omega_t'] = omega_t
        self.treebricks_dict['age_univ'] = age_univ
        self.treebricks_dict['nb_of_galaxies'] = nsub[0]
        self.treebricks_dict['nb_of_subgals'] = nsub[1]
        nmax = nsub[0]+nsub[1]
        self.treebricks_dict['nmax'] = nmax
        #initialize empty list to hold all halos
        self.treebricks_dict['galaxies'] = []

        print('nbodies:',nbodies,
              'mpart:', mpart,
              'aexp:', aexp,
              'omega_t:', omega_t,
              'age:', age_univ,
              'nsub:',nsub,
              'nmax:', nmax)

        #define the length of the box in Mpc for New Horizon
        lbox_NH=2.0*0.07*(142.8571428)*aexp
        lbox_HAGN = (142.8571428)*aexp
        halfbox=lbox_NH/2.0
        self.treebricks_dict['lbox_NH'] = lbox_NH
        self.treebricks_dict['lbox_HAGN'] = lbox_HAGN

        #def read_halo():
        t1 = time.time()
        for i in range(nmax):
            gal_dict = {}

            npart = f.read_record('i')
            gal_dict['npart'] = npart.tolist()

            #members = f.read_record('i').reshape(npart, order='F')
            members = f.read_record('i')
            gal_dict['members'] = members.tolist()

            my_number = f.read_record('i')
            gal_dict['my_number'] = my_number.tolist()

            my_timestep = f.read_record('i')
            gal_dict['my_timestep'] = my_timestep.tolist()

            level_ids = f.read_record('i')
            level_ids = level_ids.tolist()
            #print('level ids',level_ids)
            level,host_gal,host_subgal,nchild,nextsub = level_ids[0],level_ids[1],level_ids[2],level_ids[3],level_ids[4]
            gal_dict['level'] = level
            gal_dict['host_gal']  = host_gal
            gal_dict['host_subgal'] = host_subgal
            gal_dict['nchild'] = nchild
            gal_dict['nextsub'] = nextsub

            mass = f.read_record('d') #d for NH
            #print('mass',mass)
            gal_dict['mass'] = mass.tolist()

            p = f.read_record('d')
            #print('p',p)
            p = p.tolist()
            py,px,pz = p[0],p[1],p[2]
            gal_dict['px'] = px
            gal_dict['py'] = py
            gal_dict['pz'] = pz

            v = f.read_record('d')
            #print('v',v)
            v = v.tolist()
            vx,vy,vz = v[0],v[1],v[2]
            gal_dict['vx'] = vx
            gal_dict['vy'] = vy
            gal_dict['vz'] = vz

            L = f.read_record('d')
            #print('L',L)
            L = L.tolist()
            Lx,Ly,Lz = L[0],L[1],L[2]
            gal_dict['Lx'] = Lx
            gal_dict['Ly'] = Ly
            gal_dict['Lz'] = Lz

            shape = f.read_record('d')
            #print('shape',shape)
            shape = shape.tolist()
            rmax,a,b,c = shape[0],shape[1],shape[2],shape[3]
            gal_dict['rmax'] = rmax
            gal_dict['a'] = a
            gal_dict['b'] = b
            gal_dict['c'] = c

            energy = f.read_record('d')
            #print('energy',energy)
            energy = energy.tolist()
            ek,ep,et = energy[0],energy[1],energy[2]
            gal_dict['ek'] = ek
            gal_dict['ep'] = ep
            gal_dict['et'] = et

            spin = f.read_record('d')
            #print('spin',spin)
            gal_dict['spin'] = spin.tolist()

            sigma = f.read_record('d')
            #print('sigma',sigma)
            sigma = sigma.tolist()
            sig,sigma_bulge,m_bulge = sigma[0],sigma[1],sigma[2]
            gal_dict['sigma'] = sig
            gal_dict['sigma_bulge'] = sigma_bulge
            gal_dict['m_bulge'] = m_bulge

            virial = f.read_record('d')
            virial = virial.tolist()
            rvir,mvir,tvir,cvel = virial[0],virial[1],virial[2],virial[3]
            gal_dict['rvir'] = rvir
            gal_dict['mvir'] = mvir
            gal_dict['tvir'] = tvir
            gal_dict['cvel'] = cvel

            halo_profile = f.read_record('d')
            #print('halo profile',halo_profile)
            halo_profile = halo_profile.tolist()
            rho_0, r_c = halo_profile[0],halo_profile[1]
            gal_dict['rho_0'] = rho_0
            gal_dict['r_c'] = r_c

            #stellar density profiles

            nbins = f.read_record('i')
            #print('nbins:',nbins)

            rr = f.read_record('d') #array (nbins)
            gal_dict['rr'] = rr
            rho = f.read_record('d')
            gal_dict['rr'] = rho

            #print('shape of density profile',np.shape(rho))

            #add gal dict to main dict

            #return halo_dict
            self.treebricks_dict['galaxies'].append(gal_dict)

        #t1 = time.time()
        #for i in range(nhaloes):
            #halo_dict = read_halo()
            #self.treebricks_dict['haloes'].append(halo_dict)

        t2 = time.time()
        print('Reading galaxies took {:0.2f} secs.'.format(t2-t1))
        print('Total time was {:0.2f} secs.'.format(t2-t0))

        return self.treebricks_dict

class ReadGalaxy:
    def __init__(self,file_path=None):

        """
        Make a galaxies dictionary out of binary file
        """
        self.file_path = file_path
        self.galaxies_dict = None
        self.read_data()

    def read_data(self):
        t0 = time.time()
        f = FortranFile(self.file_path, 'r')

        self.galaxies_dict = {}

        my_number = f.read_record('i')
        self.galaxies_dict['my_number'] = my_number

        #gal or sub-gal => merger detection?
        level = f.read_record('i')
        self.galaxies_dict['level'] = level

        #galactic mass
        mass = f.read_record('d')
        self.galaxies_dict['mass'] = mass

        #galactic position
        p = f.read_record('d')
        p = p.tolist()
        py,px,pz = p[0],p[1],p[2]
        self.galaxies_dict['px'] = px
        self.galaxies_dict['py'] = py
        self.galaxies_dict['pz'] = pz

        #galactic velocity
        v = f.read_record('d')
        v = v.tolist()
        vx,vy,vz = v[0],v[1],v[2]
        self.galaxies_dict['vx'] = vx
        self.galaxies_dict['vy'] = vy
        self.galaxies_dict['vz'] = vz

        #galactic ang mom
        L = f.read_record('d')
        L = L.tolist()
        Lx,Ly,Lz = L[0],L[1],L[2]
        self.galaxies_dict['Lx'] = Lx
        self.galaxies_dict['Ly'] = Ly
        self.galaxies_dict['Lz'] = Lz

        #number of stars
        nstars = f.read_record('i')
        self.galaxies_dict['nstars'] = nstars

        #stellar positions
        x_stars = f.read_record('d')
        self.galaxies_dict['x_stars'] = x_stars.tolist()

        y_stars = f.read_record('d')
        self.galaxies_dict['y_stars'] = y_stars.tolist()

        z_stars = f.read_record('d')
        self.galaxies_dict['z_stars'] = z_stars.tolist()

        #stellar velocity
        vx_stars = f.read_record('d')
        self.galaxies_dict['vx_stars'] = vx_stars.tolist()

        vy_stars = f.read_record('d')
        self.galaxies_dict['vy_stars'] = vy_stars.tolist()

        vz_stars = f.read_record('d')
        self.galaxies_dict['vz_stars'] = vz_stars.tolist()

        #stellar mass
        mass_stars = f.read_record('d')
        self.galaxies_dict['mass_stars'] = mass_stars.tolist()

        #star id
        id_stars = f.read_record('i')
        self.galaxies_dict['id_stars'] = id_stars.tolist()

        #stellar age
        age_stars = f.read_record('d')
        self.galaxies_dict['age_stars'] = age_stars.tolist()

        #metallicity
        zz_stars = f.read_record('d')
        self.galaxies_dict['zz_stars'] = zz_stars.tolist()

        t1 = time.time()

        print('Reading galaxy took {:0.2f} secs.'.format(t1-t0))
        return self.galaxies_dict

class ReadDat:
    def __init__(self,file_path=None):

        """
        Return a cube array
        """
        self.file_path = file_path
        self.cube = None
        self.read_cube_data()

    def read_cube_data(self):
        t0 = time.time()
        f = FortranFile(self.file_path, 'r')
        n = f.read_record('i')
        nx,ny,nz = n[0],n[1],n[2]

        self.cube = f.read_record('f4').reshape(nx,ny,nz, order='F')
        t1 = time.time()

        print('Reading cube took {:0.2f} secs.'.format(t1-t0))
        return self.cube

class ReadFilament:
    def __init__(self,file_path=None):

        """
        Make a filaments dictionary out of ASCII NDSKL file
        """
        self.file_path = file_path
        self.filament_dict = None
        self.read_data()

    def read_data(self):
        t0 = time.time()
        self.filament_dict = {}

        #read the file first and write each line into data
        data = []
        f = open(self.file_path,'r')
        for line in f:
            data.append(line)
        f.close()

        def convert_to_list(ascii_chars,type=float):
            #strip whitespace from ends
            ascii_chars = str(ascii_chars)

            char_list = list(ascii_chars.split(" "))
            char_list = ' '.join(char_list).split()

            py_list = list(map(type, char_list))

            return py_list

        header1 = data[0]
        print('header1,',header1)

        ndims = data[1]
        print('ndims,', ndims)

        comments = data[2]
        print('Comments,',comments)

        extent = data[3]
        print('Bounding box,', extent)

        #data[4] is str(Critical Points)

        ncrit = int(data[5])
        print('ncrit,', ncrit)
        self.filament_dict['ncrit'] = ncrit

        #store all data for critical points in here
        self.filament_dict['critical_points'] = []

        ##### CPs

        add_to_idx = 6
        for i in range(ncrit):
            cp_dict = {}
            i = 0

            i += add_to_idx #make sure you are at the right line in the data list
            critical_vals = data[i]

            c_idx, px, py, pz, value, pairID, boundary = convert_to_list(critical_vals)
            #next line in data

            cp_dict['cp_idx']  = c_idx
            cp_dict['px'] = px
            cp_dict['py'] = py
            cp_dict['pz'] = pz
            cp_dict['pair_ID'] = pairID
            cp_dict['boundary'] = boundary

            i += 1
            nfil = int(data[i])
            cp_dict['nfil'] = nfil
            cp_dict['destID,filID'] = []

            for k in range(nfil):

                i += 1
                cp_on_fil = data[i]
                destID, filID = convert_to_list(cp_on_fil,int)

                cp_dict['destID,filID'].append([destID,filID])

            #make this to fill out later
            cp_dict['Field Vals'] = []
            #add all info to cp dict
            self.filament_dict['critical_points'].append(cp_dict)

            add_to_idx = i + 1

        ##### Filaments

        fil_idx = i + 1
        nfils = int(data[fil_idx+1])
        self.filament_dict['nfils'] = nfils
        print('nfils,', nfils)

        #store all data for filaments in here
        self.filament_dict['filaments'] = []

        fil_add = fil_idx+2
        for i in range(nfils):
            i = 0
            fil_dict = {}

            i += fil_add #make sure you are at the right line in the data list
            fil_info = data[i]

            cp1_idx, cp2_idx, nsamp = convert_to_list(fil_info)
            nsamp = int(nsamp)

            fil_dict['cp1_idx'] = cp1_idx
            fil_dict['cp2_idx'] = cp2_idx
            fil_dict['nsamp'] = nsamp
            fil_dict['px,py,pz'] = []

            for k in range(nsamp):

                i += 1
                positions = data[i]
                px,py,pz = convert_to_list(positions)
                #print('px,py,pz:',px,py,pz)
                fil_dict['px,py,pz'].append([px,py,pz])

            #make this to fill out later
            fil_dict['Field Vals'] = []

            #add filament info to dict
            self.filament_dict['filaments'].append(fil_dict)
            fil_add = i + 1

        cp_dat_idx = i + 1

        #Field Data
        print('Reading data fields:')
        nb_cp_dat_fields = int(data[cp_dat_idx+1])
        cp_dat_add = cp_dat_idx+2
        self.filament_dict['nb_CP_fields'] = nb_cp_dat_fields
        self.filament_dict['CP_fields'] = []

        for i in range(nb_cp_dat_fields):
            i = 0
            i += cp_dat_add #make sure you are at the right line in the data list
            cp_field_info = data[i]
            print('CP field:',cp_field_info)
            self.filament_dict['CP_fields'].append(cp_field_info)

            cp_dat_add = i + 1

        cp_field_val_idx = i + 1

        cp_val_add = cp_field_val_idx
        cp_field_vals = []
        for i in range(ncrit):
            i = 0
            i += cp_val_add #make sure you are at the right line in the data list
            cp_field_val_info = data[i]
            list_of_cp_vals = convert_to_list(cp_field_val_info)
            cp_field_vals.append(list_of_cp_vals)

            cp_val_add = i + 1

        fil_dat_idx = i + 1

        #put the field vals in the right place in the dictionary
        for j in range(ncrit):
            self.filament_dict['critical_points'][j]['Field Vals'] = cp_field_vals[j]

        nb_fil_dat_fields = int(data[fil_dat_idx+1])
        self.filament_dict['nb_fil_fields'] = nb_fil_dat_fields
        self.filament_dict['fil_fields'] = []
        fil_dat_add = fil_dat_idx+2
        for i in range(nb_fil_dat_fields):
            i = 0
            i += fil_dat_add #make sure you are at the right line in the data list
            fil_field_info = data[i]
            print('Filament field:',fil_field_info)
            self.filament_dict['fil_fields'].append(fil_field_info)

            fil_dat_add = i + 1

        fil_field_val_idx = i + 1

        fil_val_add = fil_field_val_idx

        fil_field_vals = []
        for i in range(nfils):
            i = 0
            i += fil_val_add #make sure you are at the right line in the data list
            fil_field_val_info = data[i]
            list_of_fil_vals = convert_to_list(fil_field_val_info)
            fil_field_vals.append(list_of_fil_vals)

            fil_val_add = i + 1
        #put the field vals in the right place in the dictionary
        for j in range(nfils):
            self.filament_dict['filaments'][j]['Field Vals'] = fil_field_vals[j]

        t1 = time.time()
        print('Reading filaments took {:0.2f} secs.'.format(t1 - t0))

        return self.filament_dict

reason: import numpy,from scipy


text:
#!/usr/bin/env python
import numpy as np
import sys

def outliers_modified_z_score(ys):
    threshold = 4.5

    median_y = np.median(ys)
    median_absolute_deviation_y = np.median([np.abs(y - median_y) for y in ys])
    modified_z_scores = [0.6745 * (y - median_y) / median_absolute_deviation_y for y in ys]
    return np.where(np.abs(modified_z_scores) > threshold)

def main(argv=sys.argv[1:]):
    try:
        with open(argv[0]) as f:
            lines = f.read().splitlines()
            #print ("data\n{}".format(lines))
            outliers = outliers_modified_z_score(ys=np.asarray(lines, dtype=np.float32))
            print ("the outliers we have in file {0}".format(argv[0]))
            print(" outliers:\n{0}".format(outliers))
    except Exception as e:
        print ("fail to find outliers {}".format(e))

if __name__ == '__main__':
    sys.exit(main(sys.argv[1:]))


text:
import multiprocessing as mp
import time
import datetime
import serial
from google.cloud import bigtable
from multiprocessing import Process, Queue
from google.cloud import bigtable
from google.cloud.bigtable import column_family
from google.cloud.bigtable import row_filters
import numpy as np
import matplotlib.pyplot as plt

def ard_tuple(tuple_queue):
	print("Getting data from Arduino")
	ser = serial.Serial("/dev/cu.usbmodem14201", 9600)
	for i in range(500):
		data = ser.readline()
		data = data.decode("utf-8").strip("\n")
		temp_tuple = data.split("\t")
		temp_tuple = list(temp_tuple)
		timestamp = datetime.datetime.utcnow()
		temp_tuple.append(timestamp)
		temp_tuple = tuple(temp_tuple)
		tuple_queue.put(temp_tuple)

	ser.close()
	print("All data pushed to tuple_queue, exiting ard_read process")

def tuple_BT(tuple_queue):

	# configure the connection to bigtable
	print("pushing to BT")
	client = bigtable.Client(project='sensorray', admin=True)
	instance = client.instance('instance')
	table = instance.table('table')
	for _ in range(2):
		for i in range(10):
			for j in range(25):
				podID = str(i).zfill(2) + str(j).zfill(2)
				id = podID
				# read the data from queue
				temp,light,hum,moisture,timestamp = tuple_queue.get()

				data = {
					'temp': temp,
					'light' :light,
					'humidity': hum,
					'moisture'  : moisture
				}

				# write the data
				rows = []

				for key, value in data.items():
					row = table.row(id)
					row.set_cell('sensor',key,int(value),timestamp)
					rows.append(row)

				table.mutate_rows(rows)
				print("pushed "+podID+" to BT")

def get_data(sensortype):
    client = bigtable.Client(project='sensorray', admin=True)
    instance = client.instance('instance')
    table = instance.table('table')

    row_filter = row_filters.CellsColumnLimitFilter(1)
    print("Getting 250 most recent records for "+sensortype)
    slist = []
    for i in range(0,10):
        for j in range(0,25):
            pod = str(i).zfill(2) + str(j).zfill(2)
            key = pod.encode()
            row = table.read_row(key, row_filter)
            cell = row.cells['sensor'][sensortype.encode()][0]
            slist.append(int.from_bytes(cell.value, 'big'))
    slist = np.array(slist).reshape(10,25)
    return slist

if __name__ == "__main__":
	# used between two processes ard to tuple and tuple to BT
	tuple_queue = Queue()
	ar_process = mp.Process(target=ard_tuple,args=(tuple_queue,))
	ar_process.start()

	time.sleep(1)
	bt_process = mp.Process(target = tuple_BT,args=(tuple_queue,))
	bt_process.start()

	ar_process.join()
	bt_process.join()
	print("All data pushed to cloud")


text:
import tensorflow as tf
import numpy as np
import argparse

import _init_paths
from utils import Avg, tf_acc, to1hot, cross_entropy
from networks import Baseline_CNN
from dataset import MiniImageNet, Cifar, MNIST, SubDataSet

def parse_args():
    parser = argparse.ArgumentParser(description='baseline method')
    parser.add_argument('--d', dest='dataset', default='cifar')
    parser.add_argument('--m', dest='model_name', default='b84_cifar')
    parser.add_argument('--ss', dest='initial_step', default=0)
    parser.add_argument('--mi', dest='max_iter', default=None)
    parser.add_argument('--me', dest='max_epoch', type=int, default=100)
    parser.add_argument('--bs', dest='batch_size', default=64)
    parser.add_argument('--lr', dest='initial_lr', default=1e-3)
    parser.add_argument('--fm', dest='from_ckpt', default=False)
    parser.add_argument('--nw', dest='nway', type=int, default=5)
    parser.add_argument('--ks', dest='kshot', type=int, default=1)
    args = parser.parse_args()
    return args

if __name__=='__main__':
    args = parse_args()
    print ('-'*50)
    print ('args::')
    # data size : 64 * 500
    args.max_iter = int(args.max_epoch * (64 * 500 / args.batch_size))
    for arg in vars(args):
        print ('%15s : %s'%(arg, getattr(args, arg)))
    print ('-'*50)

    x_ph = tf.placeholder(tf.float32, [None,None,None,3])
    y_ph = tf.placeholder(tf.float32, [None,100])
    lr_ph = tf.placeholder(tf.float32)
    isTr = tf.constant(True)

    model = Baseline_CNN(args.model_name)
    feature = model.network(x_ph, isTr)
    logits = model.get_logit(feature)
    y = tf.nn.softmax(logits)

    loss = cross_entropy(y_ph, y)
    opt = tf.train.AdamOptimizer(lr_ph)
    train_op = opt.minimize(loss)
    update_op = tf.get_collection(tf.GraphKeys.UPDATE_OPS)
    train_step = tf.group(train_op, update_op)
    acc = tf_acc(y_ph, y)

    logit2 = model.get_another_logit(feature)
    y2 = tf.nn.softmax(logit2)
    loss2 = cross_entropy(y_ph, y2)
    opt2 = tf.train.AdamOptimizer(lr_ph, name='opt2')
    train_op2 = opt2.minimize(loss2)
    train_step2 = tf.group(train_op, update_op)
    acc2 = tf_acc(y_ph, y2)

    sess = tf.Session()
    sess.run(tf.global_variables_initializer())

    avg = Avg(desc=['train acc', 'train loss', 'lr', 'test acc'])

    if args.dataset=='imgnet':
        train_data = MiniImageNet('train')
        val_data = MiniImageNet('train_val')
    if args.dataset=='cifar':
        train_data = Cifar('cifar100', 'train')
        val_data = Cifar('cifar100', 'test')
    if args.dataset=='mnist':
        train_data = MNIST('notMNIST', 'train')
        val_data = MNIST('notMNIST', 'test')

    model.saver_init()
    if args.from_ckpt:
        model.saver_load(sess)

    for i in range(1+args.initial_step, 1+args.max_iter):
        batch_x, batch_y = train_data.next_batch(args.batch_size)
        lr = (1e-0 if i < args.max_iter * 2.0 / 3.0 else 1e-1) * args.initial_lr
        fd = {x_ph: batch_x, y_ph: batch_y, lr_ph: lr}

        if sess.run(isTr):
            p1, p2, _ = sess.run([acc, loss, train_step], fd)
        else:
            p1, p2 = sess.run([acc, loss], fd)
        avg.add(p1, 0)
        avg.add(p2, 1)
        if i % 100 == 0:
            batch_x, batch_y = val_data.next_batch(128)
            fd = {x_ph: batch_x, y_ph: batch_y}
            avg.add(1.0/lr, 2)
            avg.add(sess.run(acc, fd), 3)
            avg.show(i)

        if i % 5000 == 0:
            avg.description()
            model.saver_save(sess)

    # fine tune the model
#    meta_test_Dtrain = MiniImageNet('test')
#    meta_test_Dtest = MiniImageNet('test_val')
#
##    random_class_nway = np.random.choice(20,
##            size=args.nway, replace=False)
#    random_class_nway = [0,6,5,14,11]
#    Dtrain = SubDataSet(random_class_nway,
#            (meta_test_Dtrain.x, meta_test_Dtrain.y), k=1000)
#    Dtest = SubDataSet(random_class_nway,
#            (meta_test_Dtest.x, meta_test_Dtest.y))
#    for i in range(1, 5001):
#        lr = 1e-4
#        batch_x, batch_y = Dtrain.next_batch(32)
#        fd = {x_ph: batch_x, y_ph:batch_y, lr_ph: lr}
#        p1, p2, _  = sess.run([acc2, loss2, train_step2], fd)
#        avg.add(p1, 0)
#        avg.add(p2, 1)
#        if i % 100 == 0:
#            batch_x, batch_y = Dtest.next_batch(100)
#            fd = {x_ph: batch_x, y_ph:batch_y, lr_ph: lr}
#            p1 = sess.run(acc2, fd)
#            avg.add(p1, 2)
#            avg.show(i)
#
##    model.saver_load(sess)


text:
import tkinter
from tkinter import *
import cv2
import PIL.Image, PIL.ImageTk
import time
import argparse
import os
from keras import backend as K
from scipy.spatial import distance as dist
from imutils.video import VideoStream
from imutils import face_utils
from threading import Thread
import numpy as np
import playsound
import argparse
import imutils
import time
import dlib
import cv2
import pyttsx3
from tensorflow.keras.applications.mobilenet_v2 import preprocess_input
from collections import deque

import os
import sys
import cv2
import numpy as np
from PIL import Image
from matplotlib import pyplot as plt
import tensorflow.compat.v1 as tf
tf.disable_v2_behavior()

from tensorflow.keras.models import load_model

class TrafficSignService:

	def __init__(self, args):

		# build models

		print("[INFO] loading traffic sign detection model...")
		self.graph = tf.Graph()
		with graph.as_default():
			graph_def = tf.GraphDef()
			with tf.gfile.GFile("models/regulation_detector_frozen.pb", 'rb') as file:
				graph_def.ParseFromString(file.read())
				tf.import_graph_def(graph_def, name="")

		print("[INFO] loading traffic sign classification model...")
		self.classifier = load_model('models/traffic_sign_classifier.hdf5')

		# variables

		self.traffic_sign_relevant_states = {
			0: 'speed_20',
			1: 'speed_30',
			2: 'speed_50',
			3: 'speed_60',
			4: 'speed_70',
			5: 'speed_80',
			6: 'speed_90',
			14: 'stop',
			27: 'crosswalk'
		}

	def process(self, frame):

		# detect

		with self.graph.as_default():
			with tf.Session(graph=self.graph) as sess:
				image_tensor = self.graph.get_tensor_by_name('image_tensor:0')
				detection_boxes = self.graph.get_tensor_by_name('detection_boxes:0')
				detection_scores = self.graph.get_tensor_by_name('detection_scores:0')
				detection_classes = self.graph.get_tensor_by_name('detection_classes:0')
				num_detections = self.graph.get_tensor_by_name('num_detections:0')

				image_exp = frame[None, ...]
				im_height, im_width = frame.shape[:2]

				boxes, scores, classes, num = sess.run(
				  [detection_boxes, detection_scores, detection_classes, num_detections],
				  feed_dict={image_tensor: image_exp})

				crops = []

				for i in range(len(boxes[0])):
					xmin, ymin = int(boxes[0,0,1] * im_width), int(boxes[0,0,0] * im_height)
					xmax, ymax = int(boxes[0,0,3] * im_width), int(boxes[0,0,2] * im_height)
					cx, cy = (xmin + xmax) // 2, (ymin + ymax) // 2

					crop = sess.run(tf.image.crop_to_bounding_box(frame, ymin, xmin, ymax - ymin, xmax - xmin))
					crop = self.preprocess_img(crop)
					crops.append(crop)

		# classify

		traffic_signs = []

		for img in crops:
			img = cv2.resize(img, (48, 48), interpolation=cv2.INTER_AREA)
			pred = self.classifier.predict(np.expand_dims(img, 0))[0]
			if pred in self.traffic_sign_relevant_states:
				traffic_sign = self.traffic_sign_relevant_states[pred]
				print(f"Warning! Traffic sign '{traffic_sign}' detected.")
				traffic_signs.append(traffic_sign)

		return traffic_signs

max_stars_repo_name: MikhailKitikov/DrivingMonitor

max_issues_repo_name: MikhailKitikov/DrivingMonitor

max_forks_repo_name: MikhailKitikov/DrivingMonitor


text:
import numpy as np
import os
import time
import random
import glob
import torch
import pathlib
import torch.optim as optim
from torch.utils.data import DataLoader
from tensorboardX import SummaryWriter
from tqdm import tqdm
from nets.load_net import gnn_model
from data.data import LoadData
from train.config import Config
from train.train_molecules_graph_regression import train_epoch, evaluate_network

class DotDict(dict):
    def __init__(self, **kwds):
        self.update(kwds)
        self.__dict__ = self

def gpu_setup(use_gpu, gpu_id):
    os.environ["CUDA_DEVICE_ORDER"] = "PCI_BUS_ID"
    os.environ["CUDA_VISIBLE_DEVICES"] = str(gpu_id)

    if torch.cuda.is_available() and use_gpu:
        print('cuda available with GPU:', torch.cuda.get_device_name(0))
        device = torch.device("cuda")
    else:
        print('cuda not available')
        device = torch.device("cpu")
    return device

def view_model_param(MODEL_NAME, net_params):
    model = gnn_model(MODEL_NAME, net_params)
    total_param = 0
    print("MODEL DETAILS:\n")
    # print(model)
    for param in model.parameters():
        total_param += np.prod(list(param.data.size()))
    print('MODEL/Total parameters:', MODEL_NAME, total_param)
    return total_param

def train_val_pipeline(model_name, dataset, params, net_params, dirs):
    t0 = time.time()
    per_epoch_time = []

    dataset_name = dataset.name

    if net_params['lap_pos_enc']:
        st = time.time()
        print("[!] Adding Laplacian positional encoding.")
        dataset._add_laplacian_positional_encodings(net_params['pos_enc_dim'])
        print('Time LapPE:', time.time() - st)

    if net_params['wl_pos_enc']:
        st = time.time()
        print("[!] Adding WL positional encoding.")
        dataset._add_wl_positional_encodings()
        print('Time WL PE:', time.time() - st)

    if net_params['full_graph']:
        st = time.time()
        print("[!] Converting the given graphs to full graphs..")
        dataset._make_full_graph()
        print('Time taken to convert to full graphs:', time.time() - st)

    trainset, valset, testset = dataset.train, dataset.val, dataset.test

    root_log_dir, root_ckpt_dir, write_file_name, write_config_file = dirs
    device = net_params['device']

    # Write the network and optimization hyper-parameters in folder config/
    with open(write_config_file + '.txt', 'w') as f:
        f.write("""Dataset: {},\nModel: {}\n\nparams={}\n\nnet_params={}\n\n\nTotal Parameters: {}\n\n""".format(
            dataset_name, model_name, params, net_params, net_params['total_param']))

    log_dir = os.path.join(root_log_dir, "RUN_" + str(0))
    writer = SummaryWriter(log_dir=log_dir)

    # setting seeds
    random.seed(params['seed'])
    np.random.seed(params['seed'])
    torch.manual_seed(params['seed'])
    if device.type == 'cuda':
        torch.cuda.manual_seed(params['seed'])

    print("Training Graphs: ", len(trainset))
    print("Validation Graphs: ", len(valset))
    print("Test Graphs: ", len(testset))

    model = gnn_model(model_name, net_params)
    model = model.to(device)

    optimizer = optim.Adam(model.parameters(), lr=params['init_lr'], weight_decay=params['weight_decay'])
    scheduler = optim.lr_scheduler.ReduceLROnPlateau(optimizer, mode='min',
                                                     factor=params['lr_reduce_factor'],
                                                     patience=params['lr_schedule_patience'],
                                                     verbose=True)

    epoch_train_losses, epoch_val_losses = [], []
    epoch_train_MAEs, epoch_val_MAEs = [], []

    train_loader = DataLoader(trainset, batch_size=params['batch_size'], shuffle=True, collate_fn=dataset.collate)
    val_loader = DataLoader(valset, batch_size=params['batch_size'], shuffle=False, collate_fn=dataset.collate)
    test_loader = DataLoader(testset, batch_size=params['batch_size'], shuffle=False, collate_fn=dataset.collate)

    # At any point you can hit Ctrl + C to break out of training early.
    try:
        with tqdm(range(params['epochs'])) as t:
            for epoch in t:

                t.set_description('Epoch %d' % epoch)

                start = time.time()

                epoch_train_loss, epoch_train_mae, optimizer = train_epoch(model, optimizer, device, train_loader)

                epoch_val_loss, epoch_val_mae = evaluate_network(model, device, val_loader)
                _, epoch_test_mae = evaluate_network(model, device, test_loader)

                epoch_train_losses.append(epoch_train_loss)
                epoch_val_losses.append(epoch_val_loss)
                epoch_train_MAEs.append(epoch_train_mae)
                epoch_val_MAEs.append(epoch_val_mae)

                writer.add_scalar('train/_loss', epoch_train_loss, epoch)
                writer.add_scalar('val/_loss', epoch_val_loss, epoch)
                writer.add_scalar('train/_mae', epoch_train_mae, epoch)
                writer.add_scalar('val/_mae', epoch_val_mae, epoch)
                writer.add_scalar('test/_mae', epoch_test_mae, epoch)
                writer.add_scalar('learning_rate', optimizer.param_groups[0]['lr'], epoch)

                t.set_postfix(time=time.time() - start, lr=optimizer.param_groups[0]['lr'],
                              train_loss=epoch_train_loss, val_loss=epoch_val_loss,
                              train_MAE=epoch_train_mae, val_MAE=epoch_val_mae,
                              test_MAE=epoch_test_mae)

                per_epoch_time.append(time.time() - start)

                # Saving checkpoint
                ckpt_dir = os.path.join(root_ckpt_dir, "RUN_")
                if not os.path.exists(ckpt_dir):
                    os.makedirs(ckpt_dir)
                torch.save(model.state_dict(), '{}.pkl'.format(ckpt_dir + "/epoch_" + str(epoch)))

                files = glob.glob(ckpt_dir + '/*.pkl')
                for file in files:
                    epoch_nb = file.split('_')[-1]
                    epoch_nb = int(epoch_nb.split('.')[0])
                    if epoch_nb < epoch - 1:
                        os.remove(file)

                scheduler.step(epoch_val_loss)

                if optimizer.param_groups[0]['lr'] < params['min_lr']:
                    print("\n!! LR EQUAL TO MIN LR SET.")
                    break

                # Stop training after params['max_time'] hours
                if time.time() - t0 > params['max_time'] * 3600:
                    print("Max_time for training elapsed {:.2f} hours, so stopping".format(params['max_time']))
                    break

    except KeyboardInterrupt:
        print('Exiting from training early because of KeyboardInterrupt')

    _, test_mae = evaluate_network(model, device, test_loader)
    _, train_mae = evaluate_network(model, device, train_loader)
    print("Test MAE: {:.4f}".format(test_mae))
    print("Train MAE: {:.4f}".format(train_mae))
    print("Convergence Time (Epochs): {:.4f}".format(epoch))
    print("TOTAL TIME TAKEN: {:.4f}s".format(time.time() - t0))
    print("AVG TIME PER EPOCH: {:.4f}s".format(np.mean(per_epoch_time)))

    writer.close()

    """
        Write the results in out_dir/results folder
    """
    with open(write_file_name + '.txt', 'w') as f:
        f.write("""Dataset: {},\nModel: {}\n\nparams={}\n\nnet_params={}\n\n{}\n\nTotal Parameters: {}\n\n
    FINAL RESULTS\nTEST MAE: {:.4f}\nTRAIN MAE: {:.4f}\n\n
    Convergence Time (Epochs): {:.4f}\nTotal Time Taken: {:.4f} hrs\nAverage Time Per Epoch: {:.4f} s\n\n\n""" \
                .format(dataset_name, model_name, params, net_params, model, net_params['total_param'],
                        test_mae, train_mae, epoch, (time.time() - t0) / 3600, np.mean(per_epoch_time)))

def main():
    cfg = Config()
    device = gpu_setup(cfg.gpu["use"], cfg.gpu["id"])

    # model, dataset, out_dir
    chr = 21
    model_name = cfg.model
    dataset_name = cfg.dataset
    dataset = LoadData(dataset_name, cfg, chr)

    params = cfg.params
    net_params = cfg.net_params
    net_params['device'] = device
    net_params['gpu_id'] = cfg.gpu['id']
    net_params['batch_size'] = params['batch_size']

    root_log_dir = cfg.output_directory + 'logs/' + model_name + time.strftime('%Hh%Mm%Ss_on_%b_%d_%Y')
    root_ckpt_dir = cfg.output_directory + 'checkpoints/' + model_name + time.strftime('%Hh%Mm%Ss_on_%b_%d_%Y')
    write_file_name = cfg.output_directory + 'results/result_' + model_name + time.strftime('%Hh%Mm%Ss_on_%b_%d_%Y')
    write_config_file = cfg.output_directory + 'configs/config_' + model_name + time.strftime('%Hh%Mm%Ss_on_%b_%d_%Y')
    dirs = root_log_dir, root_ckpt_dir, write_file_name, write_config_file

    for dir in dirs:
        directory = os.path.dirname(cfg.output_directory + dir)
        if not os.path.exists(directory):
            print("Creating directory %s" % dir)
            pathlib.Path(dir).mkdir(parents=True, exist_ok=True)
        else:
            print("Directory %s exists" % dir)

    net_params['total_param'] = view_model_param(model_name, net_params)
    train_val_pipeline(model_name, dataset, params, net_params, dirs)

main()

max_stars_repo_name: kevinbdsouza/GraphTransHiC

max_issues_repo_name: kevinbdsouza/GraphTransHiC

max_forks_repo_name: kevinbdsouza/GraphTransHiC


text:
#%%
"""
Created on August 25 2021
Displaced Diffusion and implied volatilities

This code is purely educational and comes from "Financial Engineering" course by L.A. Grzelak
The course is based on the book “Mathematical Modeling and Computation
in Finance: With Exercises and Python and MATLAB Computer Codes”,
by C.W. Oosterlee and L.A. Grzelak, World Scientific Publishing Europe Ltd, 2019.
@author: Lech A. Grzelak
"""
import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as st
import scipy.optimize as optimize
import enum

# set i= imaginary number
i   = np.complex(0.0,1.0)

# time-step needed for differentiation
dt = 0.0001

# This class defines puts and calls
class OptionType(enum.Enum):
    CALL = 1.0
    PUT = -1.0

# Black-Scholes Call option price
def BS_Call_Option_Price(CP,S_0,K,sigma,tau,r):
    if K is list:
        K = np.array(K).reshape([len(K),1])
    d1    = (np.log(S_0 / K) + (r + 0.5 * np.power(sigma,2.0)) * tau) / (sigma * np.sqrt(tau))
    d2    = d1 - sigma * np.sqrt(tau)
    if CP == OptionType.CALL:
        value = st.norm.cdf(d1) * S_0 - st.norm.cdf(d2) * K * np.exp(-r * tau)
    elif CP == OptionType.PUT:
        value = st.norm.cdf(-d2) * K * np.exp(-r * tau) - st.norm.cdf(-d1)*S_0
    return value

# Implied volatility method
def ImpliedVolatilityBlack76_xxx(CP,marketPrice,K,T,S_0):
    func = lambda sigma: np.power(BS_Call_Option_Price(CP,S_0,K,sigma,T,0.0) - marketPrice, 1.0)
    impliedVol = optimize.newton(func, 0.2, tol=1e-9)
    #impliedVol = optimize.brent(func, brack= (0.05, 2))
    return impliedVol

# Implied volatility method
def ImpliedVolatilityBlack76(CP,marketPrice,K,T,S_0):
    # To determine initial volatility we interpolate define a grid for sigma
    # and interpolate on the inverse
    sigmaGrid = np.linspace(0.0,2.0,5000)
    optPriceGrid = BS_Call_Option_Price(CP,S_0,K,sigmaGrid,T,0.0)
    sigmaInitial = np.interp(marketPrice,optPriceGrid,sigmaGrid)
    print("Initial volatility = {0}".format(sigmaInitial))

    # Use determined input for the local-search (final tuning)
    func = lambda sigma: np.power(BS_Call_Option_Price(CP,S_0,K,sigma,T,0.0) - marketPrice, 1.0)
    impliedVol = optimize.newton(func, sigmaInitial, tol=1e-5)
    print("Final volatility = {0}".format(impliedVol))
    return impliedVol

def DisplacedDiffusionModel_CallPrice(K,P0T,beta,sigma,frwd,T):
     d1    = (np.log(frwd / (beta*K+(1.0-beta)*frwd)) + 0.5 * np.power(sigma*beta,2.0) * T) / (sigma * beta* np.sqrt(T))
     d2    = d1 - sigma * beta * np.sqrt(T)
     return P0T(T) * (frwd/beta * st.norm.cdf(d1) - (K + (1.0-beta)/beta*frwd) * st.norm.cdf(d2))

def mainCalculation():
    CP  = OptionType.CALL

    K = np.linspace(0.3,2.8,22)
    K = np.array(K).reshape([len(K),1])

    # We define a ZCB curve (obtained from the market)
    P0T = lambda T: np.exp(-0.05*T)

    # DD model parameters
    beta  = 0.5
    sigma = 0.15

    # Forward rate
    frwd = 1.0

    # Maturity
    T = 2.0

    # Effect of sigma
    plt.figure(1)
    plt.grid()
    plt.xlabel('strike, K')
    plt.ylabel('implied volatility')
    sigmaV = [0.1,0.2,0.3,0.4]
    legend = []
    for sigmaTemp in sigmaV:
       callPrice = DisplacedDiffusionModel_CallPrice(K,P0T,beta,sigmaTemp,frwd,T)
       # Implied volatilities
       IV =np.zeros([len(K),1])
       for idx in range(0,len(K)):
           valCOSFrwd = callPrice[idx]/P0T(T)
           IV[idx] = ImpliedVolatilityBlack76(CP,valCOSFrwd,K[idx],T,frwd)
       plt.plot(K,IV*100.0)
       legend.append('sigma={0}'.format(sigmaTemp))
       plt.ylim([0.0,60])
    plt.legend(legend)

    # Effect of beta
    plt.figure(2)
    plt.grid()
    plt.xlabel('strike, K')
    plt.ylabel('implied volatility')
    betaV = [-0.5, 0.0001, 0.5, 1.0]
    legend = []
    for betaTemp in betaV:
       callPrice = DisplacedDiffusionModel_CallPrice(K,P0T,betaTemp,sigma,frwd,T)
       # Implied volatilities
       IV =np.zeros([len(K),1])
       for idx in range(0,len(K)):
           valCOSFrwd = callPrice[idx]/P0T(T)
           IV[idx] = ImpliedVolatilityBlack76(CP,valCOSFrwd,K[idx],T,frwd)
       plt.plot(K,IV*100.0)
       legend.append('beta={0}'.format(betaTemp))
    plt.legend(legend)

mainCalculation()

max_stars_repo_path: Lecture 11-Market Model and Convexity Adjustments/Materials/DD_ImpliedVolatility.py

reason: import numpy,import scipy