text:
using DailyMed
using EzXML
using Test

a, meta = applicationnumbers(extra = ["page" => "132"])
@test last(a) == "VMF006011"

a, meta = drugclasses()
@test first(first(a).name) == '4'

a, meta = drugnames(extra = ["page" => "930"])
@test first(a)[1:2] == "ZI"

a, meta = ndcs(extra = ["page" => "3600"])
@test last(a) == "W4215-S1393-01"

a, meta = spls(extra = ["rxcui" => "312962"])
first(a).title == "SIM"

a, meta = spls_setid("1efe378e-fee1-4ae9-8ea5-0fe2265fe2d8")
@test contains(a, "EDECRIN")

a, meta = history("9aa7140c-012c-4ea6-866d-4732e915dab6")
first(a).spl_version == "3"

a, meta = media("1efe378e-fee1-4ae9-8ea5-0fe2265fe2d8")
first(a).mime_type == "image/jpeg"

a, meta = ndcs("1efe378e-fee1-4ae9-8ea5-0fe2265fe2d8")
a == "42571"

a, meta = packaging("l51b031a0-bc40-4159-a434-d48c3eadd2ca")
@test contains(a, "aminolevulinic acid")

a, meta = uniis(extra = ["page" => "60"])
@test length(a) > 500
@test a[end - 1].active_moiety == "ZUCCHINI"


text:
# Load Julia packages (libraries) needed for clip

using StatisticalRethinking

ProjDir = @__DIR__

# Define the Stan language model

m_5_2 = "
data {
  int N;
  vector[N] D;
  vector[N] M;
}
parameters {
  real a;
  real bM;
  real<lower=0> sigma;
}
model {
  vector[N] mu = a + bM * M;
  a ~ normal( 0 , 0.2 );
  bM ~ normal( 0 , 0.5 );
  sigma ~ exponential( 1 );
  D ~ normal( mu , sigma );
}
";

# Define the SampleModel
m5_2s = SampleModel("m5.2", m_5_2);

# Input data

m5_2_data = Dict(
  "N" => size(df, 1),
  "D" => df[:, :D],
  "M" => df[:, :M]
);

# Sample using cmdstan

rc = stan_sample(m5_2s, data=m5_2_data);

if success(rc)

  # Describe the draws

  dfa2 = read_samples(m5_2s; output_format=:dataframe)
  Particles(dfa2)

end

ext: jl


text:
# AUTO GENERATED FILE - DO NOT EDIT

export ''_fefferywatermark

"""
    ''_fefferywatermark(;kwargs...)
    ''_fefferywatermark(children::Any;kwargs...)
    ''_fefferywatermark(children_maker::Function;kwargs...)

A FefferyWaterMark component.

Keyword arguments:
- `children` (a list of or a singular dash component, string or number; optional): The content of the tab - will only be displayed if this tab is selected
- `id` (String; optional)
- `className` (String; optional)
- `content` (String; optional)
- `fontColor` (String; optional)
- `fontSize` (Real; optional)
- `gapX` (Real; optional)
- `gapY` (Real; optional)
- `loading_state` (optional): . loading_state has the following type: lists containing elements 'is_loading', 'prop_name', 'component_name'.
Those elements have the following types:
  - `is_loading` (Bool; optional): Determines if the component is loading or not
  - `prop_name` (String; optional): Holds which property is loading
  - `component_name` (String; optional): Holds the name of the component that is loading
- `rotate` (Real; optional)
- `style` (Dict; optional)
- `zIndex` (Real; optional)
"""
function ''_fefferywatermark(; kwargs...)
        available_props = Symbol[:children, :id, :className, :content, :fontColor, :fontSize, :gapX, :gapY, :loading_state, :rotate, :style, :zIndex]
        wild_props = Symbol[]
        return Component("''_fefferywatermark", "FefferyWaterMark", "feffery_utils_components", available_props, wild_props; kwargs...)
end

''_fefferywatermark(children::Any; kwargs...) = ''_fefferywatermark(;kwargs..., children = children)
''_fefferywatermark(children_maker::Function; kwargs...) = ''_fefferywatermark(children_maker(); kwargs...)

max_stars_repo_name: CNFeffery/feffery-utils-components

max_issues_repo_name: CNFeffery/feffery-utils-components

max_forks_repo_name: CNFeffery/feffery-utils-components


text:
function add_rule(line::String, line_idx::Int, rules::Array{Int,2})
    first_idx = findnext(':', line, 1) + 1
    second_idx = findnext('-', line, first_idx + 1) - 1
    rules[2 * line_idx - 1, 1] = parse(Int, SubString(line, first_idx, second_idx))
    first_idx = second_idx + 2
    second_idx = findnext(' ', line, first_idx + 1) - 1
    rules[2 * line_idx - 1, 2] = parse(Int, SubString(line, first_idx, second_idx))
    first_idx = second_idx + 5
    second_idx = findnext('-', line, first_idx + 1) - 1
    rules[2 * line_idx, 1] = parse(Int, SubString(line, first_idx, second_idx))
    first_idx = second_idx + 2
    second_idx = lastindex(line)
    rules[2 * line_idx, 2] = parse(Int, SubString(line, first_idx, second_idx))
end

function is_valid_ticket(values::Array{Int,1}, rules::Array{Int,2}, n_rules::Int)::Bool
    for v in values
        if !any(rules[i, 1] <= v <= rules[i, 2] for i in 1:n_rules)
            return false
        end
    end
    return true
end

function update_possibilities(possibilities::Array{Int,2}, rules::Array{Int,2}, values::Array{Int,1}, n_rules::Int)
    for row_idx in 1:n_rules
        for col_idx in 1:n_rules
            i = 2 * col_idx - 1
            if !(rules[i, 1] <= values[row_idx] <= rules[i, 2]) && !(rules[i + 1, 1] <= values[row_idx] <= rules[i + 1, 2])
                possibilities[row_idx, col_idx] = 0
            end
        end
    end
end

function run(s)
    rules::Array{Int,2} = zeros(Int, (50, 2))
    rules_mask::Array{Bool,1} = zeros(Bool, (0))
    n_rules::Int = 0

    myticket::Array{Int,1} = zeros(Int, (0))
    lines = eachline(IOBuffer(s))
    for line in lines
        if isempty(line)
            break
        end
        n_rules += 1
        add_rule(line, n_rules, rules)
        push!(rules_mask, startswith(line, "departure"))
    end
    possibilities::Array{Int,2} = ones(Int, (n_rules, n_rules))
    for line in Iterators.drop(lines, 1)
        if isempty(line)
            break
        end
        myticket = parse.(Int, split(line, ","))
        update_possibilities(possibilities, rules, myticket, n_rules)
    end
    for line in Iterators.drop(lines, 1)
        values = parse.(Int, split(line, ","))
        if is_valid_ticket(values, rules, n_rules)
            update_possibilities(possibilities, rules, values, n_rules)
        end
    end
    while sum(possibilities) > n_rules
        for col_idx in 1:n_rules
            v = view(possibilities, :, col_idx)
            if sum(v) == 1
                row_idx = argmax(v)
                possibilities[row_idx, 1:end] .= 0
                possibilities[row_idx, col_idx] = 1
            end
        end
    end
    indexes = [argmax(view(possibilities, :, col_idx)) for col_idx in 1:n_rules if rules_mask[col_idx]]
    return prod(myticket[idx] for idx in indexes)
end

#########################################

function main()
    res, time, memory = @timed run(ARGS[1])
    println("_duration:$(time * 1000)")
    println(res)
end

main()


text:
################################################################################
#                             Lines, LineSegments                              #
################################################################################

function draw_atomic(scene::Scene, screen::CairoScreen, primitive::Union{Lines, LineSegments})
    fields = @get_attribute(primitive, (color, linewidth, linestyle))
    linestyle = Makie.convert_attribute(linestyle, Makie.key"linestyle"())
    ctx = screen.context
    model = primitive[:model][]
    positions = primitive[1][]

    isempty(positions) && return

    # workaround for a LineSegments object created from a GLNormalMesh
    # the input argument is a view of points using faces, which results in
    # a vector of tuples of two points. we convert those to a list of points
    # so they don't trip up the rest of the pipeline
    # TODO this shouldn't be necessary anymore!
    if positions isa SubArray{<:Point3, 1, P, <:Tuple{Array{<:AbstractFace}}} where P
        positions = let
            pos = Point3f0[]
            for tup in positions
                push!(pos, tup[1])
                push!(pos, tup[2])
            end
            pos
        end
    end

    projected_positions = project_position.(Ref(scene), positions, Ref(model))

    if color isa AbstractArray{<: Number}
        color = numbers_to_colors(color, primitive)
    end

    # color is now a color or an array of colors
    # if it's an array of colors, each segment must be stroked separately

    # The linestyle can be set globally, as we do here.
    # However, there is a discrepancy between Makie
    # and Cairo when it comes to linestyles.
    # For Makie, the linestyle array is cumulative,
    # and defines the "absolute" endpoints of segments.
    # However, for Cairo, each value provides the length of
    # alternate "on" and "off" portions of the stroke.
    # Therefore, we take the diff of the given linestyle,
    # to convert the "absolute" coordinates into "relative" ones.
    if !isnothing(linestyle) && !(linewidth isa AbstractArray)
        Cairo.set_dash(ctx, diff(Float64.(linestyle)) .* linewidth)
    end
    if color isa AbstractArray || linewidth isa AbstractArray
        # stroke each segment separately, this means disjointed segments with probably
        # wonky dash patterns if segments are short

        # we can hide the gaps by setting the line cap to round
        Cairo.set_line_cap(ctx, Cairo.CAIRO_LINE_CAP_ROUND)
        draw_multi(
            primitive, ctx,
            projected_positions,
            color, linewidth,
            isnothing(linestyle) ? nothing : diff(Float64.(linestyle))
        )
    else
        # stroke the whole line at once if it has only one color
        # this allows correct linestyles and line joins as well and will be the
        # most common case
        Cairo.set_line_width(ctx, linewidth)
        Cairo.set_source_rgba(ctx, red(color), green(color), blue(color), alpha(color))
        draw_single(primitive, ctx, projected_positions)
    end
    nothing
end

function draw_single(primitive::Lines, ctx, positions)
    n = length(positions)
    @inbounds for i in 1:n
        p = positions[i]
        # only take action for non-NaNs
        if !isnan(p)
            # new line segment at beginning or if previously NaN
            if i == 1 || isnan(positions[i-1])
                Cairo.move_to(ctx, p...)
            else
                Cairo.line_to(ctx, p...)
                # complete line segment at end or if next point is NaN
                if i == n || isnan(positions[i+1])
                    Cairo.stroke(ctx)
                end
            end
        end
    end
end

function draw_single(primitive::LineSegments, ctx, positions)

    @assert iseven(length(positions))

    @inbounds for i in 1:2:length(positions)-1
        p1 = positions[i]
        p2 = positions[i+1]

        if isnan(p1) || isnan(p2)
            continue
        else
            Cairo.move_to(ctx, p1...)
            Cairo.line_to(ctx, p2...)
            Cairo.stroke(ctx)
        end
    end
end

# if linewidth is not an array
function draw_multi(primitive, ctx, positions, colors::AbstractArray, linewidth, dash)
    draw_multi(primitive, ctx, positions, colors, [linewidth for c in colors], dash)
end

# if color is not an array
function draw_multi(primitive, ctx, positions, color, linewidths::AbstractArray, dash)
    draw_multi(primitive, ctx, positions, [color for l in linewidths], linewidths, dash)
end

function draw_multi(primitive::Union{Lines, LineSegments}, ctx, positions, colors::AbstractArray, linewidths::AbstractArray, dash)
    if primitive isa LineSegments
        @assert iseven(length(positions))
    end
    @assert length(positions) == length(colors)
    @assert length(linewidths) == length(colors)

    iterator = if primitive isa Lines
        1:length(positions)-1
    elseif primitive isa LineSegments
        1:2:length(positions)
    end

    for i in iterator
        if isnan(positions[i+1]) || isnan(positions[i])
            continue
        end
        Cairo.move_to(ctx, positions[i]...)

        Cairo.line_to(ctx, positions[i+1]...)
        if linewidths[i] != linewidths[i+1]
            error("Cairo doesn't support two different line widths ($(linewidths[i]) and $(linewidths[i+1])) at the endpoints of a line.")
        end
        Cairo.set_line_width(ctx, linewidths[i])
        !isnothing(dash) && Cairo.set_dash(ctx, dash .* linewidths[i])
        c1 = colors[i]
        c2 = colors[i+1]
        # we can avoid the more expensive gradient if the colors are the same
        # this happens if one color was given for each segment
        if c1 == c2
            Cairo.set_source_rgba(ctx, red(c1), green(c1), blue(c1), alpha(c1))
            Cairo.stroke(ctx)
        else
            pat = Cairo.pattern_create_linear(positions[i]..., positions[i+1]...)
            Cairo.pattern_add_color_stop_rgba(pat, 0, red(c1), green(c1), blue(c1), alpha(c1))
            Cairo.pattern_add_color_stop_rgba(pat, 1, red(c2), green(c2), blue(c2), alpha(c2))
            Cairo.set_source(ctx, pat)
            Cairo.stroke(ctx)
            Cairo.destroy(pat)
        end
    end
end

################################################################################
#                                   Scatter                                    #
################################################################################

function draw_atomic(scene::Scene, screen::CairoScreen, primitive::Scatter)
    fields = @get_attribute(primitive, (color, markersize, strokecolor, strokewidth, marker, marker_offset, rotations))
    @get_attribute(primitive, (transform_marker,))

    ctx = screen.context
    model = primitive[:model][]
    positions = primitive[1][]
    isempty(positions) && return
    size_model = transform_marker ? model : Mat4f0(I)

    font = to_font(to_value(get(primitive, :font, Makie.defaultfont())))

    colors = if color isa AbstractArray{<: Number}
        numbers_to_colors(color, primitive)
    else
        color
    end

    broadcast_foreach(primitive[1][], colors, markersize, strokecolor,
                      strokewidth, marker, marker_offset, remove_billboard(rotations)) do point, col,
                          markersize, strokecolor, strokewidth, marker, mo, rotation

        # if we give size in pixels, the size is always equal to that value
        is_pixelspace = haskey(primitive, :markerspace) && primitive.markerspace[] == Makie.Pixel
        scale = if is_pixelspace
            Makie.to_2d_scale(markersize)
        else
            # otherwise calculate a scaled size
            project_scale(scene, markersize, size_model)
        end
        offset = if is_pixelspace
            Makie.to_2d_scale(mo)
        else
            project_scale(scene, mo, size_model)
        end

        pos = project_position(scene, point, model)

        isnan(pos) && return

        Cairo.set_source_rgba(ctx, rgbatuple(col)...)
        m = convert_attribute(marker, key"marker"(), key"scatter"())
        Cairo.save(ctx)
        if m isa Char
            draw_marker(ctx, m, best_font(m, font), pos, scale, strokecolor, strokewidth, offset, rotation)
        else
            draw_marker(ctx, m, pos, scale, strokecolor, strokewidth, offset, rotation)
        end
        Cairo.restore(ctx)
    end
    nothing
end

function draw_marker(ctx, marker::Char, font, pos, scale, strokecolor, strokewidth, marker_offset, rotation)
    # Marker offset is meant to be relative to the
    # bottom left corner of the box centered at
    # `pos` with sides defined by `scale`, but
    # this does not take the character's dimensions
    # into account.
    # Here, we reposition the marker offset to be
    # relative to the center of the char.
    marker_offset = marker_offset .+ scale ./ 2

    cairoface = set_ft_font(ctx, font)

    charextent = Makie.FreeTypeAbstraction.internal_get_extent(font, marker)
    inkbb = Makie.FreeTypeAbstraction.inkboundingbox(charextent)

    # scale normalized bbox by font size
    inkbb_scaled = FRect2D(origin(inkbb) .* scale, widths(inkbb) .* scale)

    # flip y for the centering shift of the character because in Cairo y goes down
    centering_offset = [1, -1] .* (-origin(inkbb_scaled) .- 0.5 .* widths(inkbb_scaled))
    # this is the origin where we actually have to place the glyph so it can be centered
    charorigin = pos .+ Vec2f0(marker_offset[1], -marker_offset[2])
    old_matrix = get_font_matrix(ctx)
    set_font_matrix(ctx, scale_matrix(scale...))

    # First, we translate to the point where the
    # marker is supposed to go.
    Cairo.translate(ctx, charorigin...)
    # Then, we rotate the context by the
    # appropriate amount,
    Cairo.rotate(ctx, to_2d_rotation(rotation))
    # and apply a centering offset to account for
    # the fact that text is shown from the (relative)
    # bottom left corner.
    Cairo.translate(ctx, centering_offset...)

    Cairo.move_to(ctx, 0, 0)
    Cairo.text_path(ctx, string(marker))
    Cairo.fill_preserve(ctx)
    # stroke
    Cairo.set_line_width(ctx, strokewidth)
    Cairo.set_source_rgba(ctx, rgbatuple(strokecolor)...)
    Cairo.stroke(ctx)

    # if we use set_ft_font we should destroy the pointer it returns
    cairo_font_face_destroy(cairoface)

    set_font_matrix(ctx, old_matrix)
end

function draw_marker(ctx, marker::Circle, pos, scale, strokecolor, strokewidth, marker_offset, rotation)

    marker_offset = marker_offset + scale ./ 2
    pos += Point2f0(marker_offset[1], -marker_offset[2])
    Cairo.arc(ctx, pos[1], pos[2], scale[1]/2, 0, 2*pi)
    Cairo.fill_preserve(ctx)

    Cairo.set_line_width(ctx, Float64(strokewidth))

    sc = to_color(strokecolor)
    Cairo.set_source_rgba(ctx, rgbatuple(sc)...)
    Cairo.stroke(ctx)
end

function draw_marker(ctx, marker::Rect, pos, scale, strokecolor, strokewidth, marker_offset, rotation)
    s2 = Point2((widths(marker) .* scale .* (1, -1))...)
    pos = pos .+ Point2f0(marker_offset[1], -marker_offset[2])
    Cairo.rotate(ctx, to_2d_rotation(rotation))
    Cairo.rectangle(ctx, pos[1], pos[2], s2...)
    Cairo.fill_preserve(ctx)
    Cairo.set_line_width(ctx, Float64(strokewidth))
    sc = to_color(strokecolor)
    Cairo.set_source_rgba(ctx, rgbatuple(sc)...)
    Cairo.stroke(ctx)
end

################################################################################
#                                     Text                                     #
################################################################################

function p3_to_p2(p::Point3{T}) where T
    if p[3] == 0 || isnan(p[3])
        Point2{T}(p[1:2]...)
    else
        error("Can't reduce Point3 to Point2 with nonzero third component $(p[3]).")
    end
end

function draw_atomic(scene::Scene, screen::CairoScreen, primitive::Text{<:Tuple{<:G}}) where G <: Union{AbstractArray{<:Makie.GlyphCollection}, Makie.GlyphCollection}
    ctx = screen.context
    @get_attribute(primitive, (rotation, model, space, offset))
    position = primitive.attributes[:position][]
    # use cached glyph info
    glyph_collection = to_value(primitive[1])

    draw_glyph_collection(scene, ctx, position, glyph_collection, remove_billboard(rotation), model, space, offset)

    nothing
end

function draw_glyph_collection(scene, ctx, positions, glyph_collections::AbstractArray, rotation, model::SMatrix, space, offset)

    # TODO: why is the Ref around model necessary? doesn't broadcast_foreach handle staticarrays matrices?
    broadcast_foreach(positions, glyph_collections, rotation,
        Ref(model), space, offset) do pos, glayout, ro, mo, sp, off

        draw_glyph_collection(scene, ctx, pos, glayout, ro, mo, sp, off)
    end
end

function draw_glyph_collection(scene, ctx, position, glyph_collection, rotation, model, space, offsets)

    glyphs = glyph_collection.glyphs
    glyphoffsets = glyph_collection.origins
    fonts = glyph_collection.fonts
    rotations = glyph_collection.rotations
    scales = glyph_collection.scales
    colors = glyph_collection.colors
    strokewidths = glyph_collection.strokewidths
    strokecolors = glyph_collection.strokecolors

    Cairo.save(ctx)

    broadcast_foreach(glyphs, glyphoffsets, fonts, rotations, scales, colors, strokewidths, strokecolors, offsets) do glyph,
        glyphoffset, font, rotation, scale, color, strokewidth, strokecolor, offset

        cairoface = set_ft_font(ctx, font)
        old_matrix = get_font_matrix(ctx)

        p3_offset = to_ndim(Point3f0, offset, 0)

        glyph in ('\r', '\n') && return

        Cairo.save(ctx)
        Cairo.set_source_rgba(ctx, rgbatuple(color)...)

        if space == :data
            # in data space, the glyph offsets are just added to the string positions
            # and then projected

            # glyph position in data coordinates (offset has rotation applied already)
            gpos_data = to_ndim(Point3f0, position, 0) .+ glyphoffset .+ p3_offset

            scale3 = scale isa Number ? Point3f0(scale, scale, 0) : to_ndim(Point3f0, scale, 0)

            # this could be done better but it works at least

            # the CairoMatrix is found by transforming the right and up vector
            # of the character into screen space and then subtracting the projected
            # origin. The resulting vectors give the directions in which the character
            # needs to be stretched in order to match the 3D projection

            xvec = rotation * (scale3[1] * Point3f0(1, 0, 0))
            yvec = rotation * (scale3[2] * Point3f0(0, -1, 0))

            glyphpos = project_position(scene, gpos_data, Mat4f0(I))
            xproj = project_position(scene, gpos_data + xvec, Mat4f0(I))
            yproj = project_position(scene, gpos_data + yvec, Mat4f0(I))

            xdiff = xproj - glyphpos
            ydiff = yproj - glyphpos

            mat = Cairo.CairoMatrix(
                xdiff[1], xdiff[2],
                ydiff[1], ydiff[2],
                0, 0,
            )

        elseif space == :screen
            # in screen space, the glyph offsets are added after projecting
            # the string position into screen space
            glyphpos = project_position(
                scene,
                position,
                Mat4f0(I)) .+ (p3_to_p2(glyphoffset .+ p3_offset)) .* (1, -1) # flip for Cairo
            # and the scale is just taken as is
            scale = length(scale) == 2 ? scale : SVector(scale, scale)

            mat = scale_matrix(scale...)
        else
            error()
        end

        Cairo.save(ctx)
        Cairo.move_to(ctx, glyphpos...)
        set_font_matrix(ctx, mat)
        if space == :screen
            Cairo.rotate(ctx, to_2d_rotation(rotation))
        end
        Cairo.show_text(ctx, string(glyph))
        Cairo.restore(ctx)

        if strokewidth > 0 && strokecolor != RGBAf0(0, 0, 0, 0)
            Cairo.save(ctx)
            Cairo.move_to(ctx, glyphpos...)
            set_font_matrix(ctx, mat)
            if space == :screen
                Cairo.rotate(ctx, to_2d_rotation(rotation))
            end
            Cairo.text_path(ctx, string(glyph))
            Cairo.set_source_rgba(ctx, rgbatuple(strokecolor)...)
            Cairo.set_line_width(ctx, strokewidth)
            Cairo.stroke(ctx)
            Cairo.restore(ctx)
        end
        Cairo.restore(ctx)

        cairo_font_face_destroy(cairoface)
        set_font_matrix(ctx, old_matrix)
    end

    Cairo.restore(ctx)

    nothing
end

################################################################################
#                                Heatmap, Image                                #
################################################################################

"""
    regularly_spaced_array_to_range(arr)
If possible, converts `arr` to a range.
If not, returns array unchanged.
"""
function regularly_spaced_array_to_range(arr)
    diffs = unique!(sort!(diff(arr)))
    step = sum(diffs) ./ length(diffs)
    if all(x-> x ≈ step, diffs)
        m, M = extrema(arr)
        if step < zero(step)
            m, M = M, m
        end
        # don't use stop=M, since that may not include M
        return range(m; step=step, length=length(arr))
    else
        return arr
    end
end

regularly_spaced_array_to_range(arr::AbstractRange) = arr

"""
    interpolation_flag(is_vector, interp, wpx, hpx, w, h)

* is_vector: if we're using vector backend
* interp: does the user want to interpolate?
* wpx, hpx: projected size of the image in pixels, so the actual width in pixels on screen
* w, h: size of image in pixels
"""
function interpolation_flag(is_vector, interp, wpx, hpx, w, h)
    if interp
        if is_vector
            return Cairo.FILTER_BILINEAR
        else
            return Cairo.FILTER_BEST
        end
    else
        if wpx < w || hpx < h
            # if size of image size in pixels is larger then the rectangle it gets drawn into,
            # the pixels will be smaller than what ends up on screen, so one won't be able to see rectangles.
            # In that case, we need to apply filtering, or we get artifacts from incorrectly downsampling!
            return interpolation_flag(is_vector, true, wpx, hpx, w, h)
        else
            return Cairo.FILTER_NEAREST
        end
    end
end

function draw_atomic(scene::Scene, screen::CairoScreen, primitive::Union{Heatmap, Image})
    ctx = screen.context
    image = primitive[3][]
    xs, ys = primitive[1][], primitive[2][]
    if !(xs isa AbstractVector)
        l, r = extrema(xs)
        N = size(image, 1)
        xs = range(l, r, length = N+1)
    else
        xs = regularly_spaced_array_to_range(xs)
    end
    if !(ys isa AbstractVector)
        l, r = extrema(ys)
        N = size(image, 2)
        ys = range(l, r, length = N+1)
    else
        ys = regularly_spaced_array_to_range(ys)
    end
    model = primitive[:model][]
    interp = to_value(get(primitive, :interpolate, true))
    weird_cairo_limit = (2^15) - 23

    # Debug attribute we can set to disable fastpath
    # probably shouldn't really be part of the interface
    fast_path = to_value(get(primitive, :fast_path, true))
    # Vector backends don't support FILTER_NEAREST for interp == false, so in that case we also need to draw rects
    is_vector = is_vector_backend(ctx)
    t = Makie.transform_func_obs(primitive)[]
    identity_transform = t === identity || t isa Tuple && all(x-> x === identity, t)
    if fast_path && xs isa AbstractRange && ys isa AbstractRange && !(is_vector && !interp) && identity_transform
        imsize = ((first(xs), last(xs)), (first(ys), last(ys)))

        # find projected image corners
        # this already takes care of flipping the image to correct cairo orientation
        xy = project_position(scene, Point2f0(first.(imsize)), model)
        xymax = project_position(scene, Point2f0(last.(imsize)), model)
        w, h = xymax .- xy

        s = to_cairo_image(image, primitive)

        interp_flag = interpolation_flag(is_vector, interp, abs(w), abs(h), s.width, s.height)

        if s.width > weird_cairo_limit || s.height > weird_cairo_limit
            error("Cairo stops rendering images bigger than $(weird_cairo_limit), which is likely a bug in Cairo. Please resample your image/heatmap with e.g. `ImageTransformations.imresize`")
        end

        Cairo.rectangle(ctx, xy..., w, h)
        Cairo.save(ctx)
        Cairo.translate(ctx, xy...)
        Cairo.scale(ctx, w / s.width, h / s.height)
        Cairo.set_source_surface(ctx, s, 0, 0)
        p = Cairo.get_source(ctx)
        # this is needed to avoid blurry edges
        Cairo.pattern_set_extend(p, Cairo.EXTEND_PAD)
        # Set filter doesn't work!?
        Cairo.pattern_set_filter(p, interp_flag)
        Cairo.fill(ctx)
        Cairo.restore(ctx)

    else
        # We need to draw rectangles for vector backends, or irregular grids
        if interp
            error("Interpolation for non gridded heatmaps/images isn't supported right now. Please use interpolate=false for this plot")
        end
        # find projected image corners
        # this already takes care of flipping the image to correct cairo orientation
        xys = [project_position(scene, Point2f0(x, y), model) for x in xs, y in ys]
        colors = to_rgba_image(image, primitive)

        # Note: xs and ys should have size ni+1, nj+1
        ni, nj = size(image)
        if ni + 1 != length(xs) || nj + 1 != length(ys)
            error("Error in conversion pipeline. xs and ys should have size ni+1, nj+1. Found: xs: $(length(xs)), ys: $(length(ys)), ni: $(ni), nj: $(nj)")
        end
        @inbounds for i in 1:ni, j in 1:nj
            x0, y0 = xys[i, j]
            x1, y1 = xys[i+1, j+1]
            w = x1 - x0; h = y1 - y0

            # there are usually white lines between directly adjacent rectangles
            # in vector graphics because of anti-aliasing

            # if we let each cell stick out (bulge) a little bit (half a point) under its neighbors
            # those lines disappear

            # we heuristically only do this if the adjacent cells are fully opaque
            # and if we're not in the last row / column so the overall heatmap doesn't get bigger

            # this should be the most common case by far, though

            xbulge = if i < ni && alpha(colors[i+1, j]) == 1
                0.5
            else
                0.0
            end
            ybulge = if j < nj && alpha(colors[i, j+1]) == 1
                0.5
            else
                0.0
            end

            # we add the bulge in the direction of cell width / height in case the axes are reversed
            Cairo.rectangle(ctx, x0, y0, w + sign(w) * xbulge, h + sign(h) * ybulge)
            Cairo.set_source_rgba(ctx, rgbatuple(colors[i, j])...)
            Cairo.fill(ctx)
        end
    end
end

################################################################################
#                                     Mesh                                     #
################################################################################

function draw_atomic(scene::Scene, screen::CairoScreen, primitive::Makie.Mesh)
    if scene.camera_controls[] isa Union{Camera2D, Makie.PixelCamera}
        draw_mesh2D(scene, screen, primitive)
    else
        if !haskey(primitive, :faceculling)
            primitive[:faceculling] = Node(-10)
        end
        draw_mesh3D(scene, screen, primitive)
    end
    return nothing
end

function draw_mesh2D(scene, screen, primitive)
    @get_attribute(primitive, (color,))

    colormap = get(primitive, :colormap, nothing) |> to_value |> to_colormap
    colorrange = get(primitive, :colorrange, nothing) |> to_value
    ctx = screen.context
    model = primitive.model[]
    mesh = GeometryBasics.mesh(primitive[1][])
    # Priorize colors of the mesh if present
    # This is a hack, which needs cleaning up in the Mesh plot type!
    color = hasproperty(mesh, :color) ? mesh.color : color
    vs =  decompose(Point, mesh); fs = decompose(TriangleFace, mesh)
    uv = hasproperty(mesh, :uv) ? mesh.uv : nothing
    pattern = Cairo.CairoPatternMesh()

    cols = per_face_colors(color, colormap, colorrange, nothing, vs, fs, nothing, uv)
    for (f, (c1, c2, c3)) in zip(fs, cols)
        t1, t2, t3 =  project_position.(scene, vs[f], (model,)) #triangle points
        Cairo.mesh_pattern_begin_patch(pattern)

        Cairo.mesh_pattern_move_to(pattern, t1...)
        Cairo.mesh_pattern_line_to(pattern, t2...)
        Cairo.mesh_pattern_line_to(pattern, t3...)

        mesh_pattern_set_corner_color(pattern, 0, c1)
        mesh_pattern_set_corner_color(pattern, 1, c2)
        mesh_pattern_set_corner_color(pattern, 2, c3)

        Cairo.mesh_pattern_end_patch(pattern)
    end
    Cairo.set_source(ctx, pattern)
    Cairo.close_path(ctx)
    Cairo.paint(ctx)
    return nothing
end

function average_z(positions, face)
    vs = positions[face]
    sum(v -> v[3], vs) / length(vs)
end

nan2zero(x) = !isnan(x) * x

function draw_mesh3D(
        scene, screen, primitive;
        mesh = primitive[1][], pos = Vec4f0(0), scale = 1f0
    )
    @get_attribute(primitive, (color, shading, lightposition, ambient, diffuse,
        specular, shininess, faceculling))

    colormap = get(primitive, :colormap, nothing) |> to_value |> to_colormap
    colorrange = get(primitive, :colorrange, nothing) |> to_value
    matcap = get(primitive, :matcap, nothing) |> to_value
    # Priorize colors of the mesh if present
    color = hasproperty(mesh, :color) ? mesh.color : color

    ctx = screen.context

    model = primitive.model[]
    view = scene.camera.view[]
    projection = scene.camera.projection[]
    normalmatrix = get(
        scene.attributes, :normalmatrix, let
            i = SOneTo(3)
            transpose(inv(view[i, i] * model[i, i]))
        end
    )

    # Mesh data
    # transform to view/camera space
    func = Makie.transform_func_obs(scene)[]
    # pass func as argument to function, so that we get a function barrier
    # and have `func` be fully typed inside closure
    vs = broadcast(decompose(Point, mesh), (func,)) do v, f
        # Should v get a nan2zero?
        v = Makie.apply_transform(f, v)
        p4d = to_ndim(Vec4f0, scale .* to_ndim(Vec3f0, v, 0f0), 1f0)
        view * (model * p4d .+ to_ndim(Vec4f0, pos, 0f0))
    end
    fs = decompose(GLTriangleFace, mesh)
    uv = hasproperty(mesh, :uv) ? mesh.uv : nothing
    ns = map(n -> normalize(normalmatrix * n), normals(mesh))
    cols = per_face_colors(
        color, colormap, colorrange, matcap, vs, fs, ns, uv,
        get(primitive, :lowclip, nothing) |> to_value |> color_or_nothing,
        get(primitive, :highclip, nothing) |> to_value |> color_or_nothing,
        get(primitive, :nan_color, nothing) |> to_value |> color_or_nothing
    )

    # Liight math happens in view/camera space
    if lightposition == :eyeposition
        lightposition = scene.camera.eyeposition[]
    end
    lightpos = (view * to_ndim(Vec4f0, lightposition, 1.0))[Vec(1, 2, 3)]

    # Camera to screen space
    ts = map(vs) do v
        clip = projection * v
        @inbounds begin
            p = (clip ./ clip[4])[Vec(1, 2)]
            p_yflip = Vec2f0(p[1], -p[2])
            p_0_to_1 = (p_yflip .+ 1f0) / 2f0
        end
        p = p_0_to_1 .* scene.camera.resolution[]
        Vec3f0(p[1], p[2], clip[3])
    end

    # Approximate zorder
    zorder = sortperm(fs, by = f -> average_z(ts, f))

    # Face culling
    zorder = filter(i -> any(last.(ns[fs[i]]) .> faceculling), zorder)

    pattern = Cairo.CairoPatternMesh()
    for k in reverse(zorder)
        f = fs[k]
        t1, t2, t3 = ts[f]

        # light calculation
        c1, c2, c3 = if shading
            map(ns[f], vs[f], cols[k]) do N, v, c
                L = normalize(lightpos .- v[Vec(1,2,3)])
                diff_coeff = max(dot(L, N), 0.0)
                H = normalize(L + normalize(-v[SOneTo(3)]))
                spec_coeff = max(dot(H, N), 0.0)^shininess
                c = RGBA(c)
                new_c = (ambient .+ diff_coeff .* diffuse) .* Vec3f0(c.r, c.g, c.b) .+
                        specular * spec_coeff
                RGBA(new_c..., c.alpha)
            end
        else
            cols[k]
        end
        # debug normal coloring
        # n1, n2, n3 = Vec3f0(0.5) .+ 0.5ns[f]
        # c1 = RGB(n1...)
        # c2 = RGB(n2...)
        # c3 = RGB(n3...)

        Cairo.mesh_pattern_begin_patch(pattern)

        Cairo.mesh_pattern_move_to(pattern, t1[1], t1[2])
        Cairo.mesh_pattern_line_to(pattern, t2[1], t2[2])
        Cairo.mesh_pattern_line_to(pattern, t3[1], t3[2])

        mesh_pattern_set_corner_color(pattern, 0, c1)
        mesh_pattern_set_corner_color(pattern, 1, c2)
        mesh_pattern_set_corner_color(pattern, 2, c3)

        Cairo.mesh_pattern_end_patch(pattern)
    end
    Cairo.set_source(ctx, pattern)
    Cairo.close_path(ctx)
    Cairo.paint(ctx)
    return nothing
end

################################################################################
#                                   Surface                                    #
################################################################################

function draw_atomic(scene::Scene, screen::CairoScreen, primitive::Makie.Surface)
    # Pretend the surface plot is a mesh plot and plot that instead
    mesh = surface2mesh(primitive[1][], primitive[2][], primitive[3][])
    old = primitive[:color]
    if old[] === nothing
        primitive[:color] = primitive[3]
    end
    if !haskey(primitive, :faceculling)
        primitive[:faceculling] = Node(-10)
    end
    draw_mesh3D(scene, screen, primitive, mesh=mesh)
    primitive[:color] = old
    return nothing
end

function surface2mesh(xs::Makie.ClosedInterval, ys, zs::AbstractMatrix)
    surface2mesh(range(xs.left, xs.right, length = size(zs, 1)), ys, zs)
end

function surface2mesh(xs, ys::Makie.ClosedInterval, zs::AbstractMatrix)
    surface2mesh(xs, range(ys.left, ys.right, length = size(zs, 2)), zs)
end

function surface2mesh(xs::Makie.ClosedInterval, ys::Makie.ClosedInterval, zs::AbstractMatrix)
    surface2mesh(
        range(xs.left, xs.right, length = size(zs, 1)),
        range(ys.left, ys.right, length = size(zs, 2)),
        zs)
end

function surface2mesh(xs::AbstractVector, ys::AbstractVector, zs::AbstractMatrix)
    ps = [nan2zero.(Point3f0(xs[i], ys[j], zs[i, j])) for j in eachindex(ys) for i in eachindex(xs)]
    idxs = LinearIndices(size(zs))
    faces = [
        QuadFace(idxs[i, j], idxs[i+1, j], idxs[i+1, j+1], idxs[i, j+1])
        for j in 1:size(zs, 2)-1 for i in 1:size(zs, 1)-1
    ]
    normal_mesh(ps, faces)
end

function surface2mesh(xs::AbstractMatrix, ys::AbstractMatrix, zs::AbstractMatrix)
    ps = [nan2zero.(Point3f0(xs[i, j], ys[i, j], zs[i, j])) for j in 1:size(zs, 2) for i in 1:size(zs, 1)]
    idxs = LinearIndices(size(zs))
    faces = [
        QuadFace(idxs[i, j], idxs[i+1, j], idxs[i+1, j+1], idxs[i, j+1])
        for j in 1:size(zs, 2)-1 for i in 1:size(zs, 1)-1
    ]
    normal_mesh(ps, faces)
end

################################################################################
#                                 MeshScatter                                  #
################################################################################

function draw_atomic(scene::Scene, screen::CairoScreen, primitive::Makie.MeshScatter)
    @get_attribute(primitive, (color, model, marker, markersize, rotations))

    if color isa AbstractArray{<: Number}
        color = numbers_to_colors(color, primitive)
    end

    m = convert_attribute(marker, key"marker"(), key"meshscatter"())
    pos = primitive[1][]
    # For correct z-ordering we need to be in view/camera or screen space
    model = copy(model)
    view = scene.camera.view[]

    zorder = sortperm(pos, by = p -> begin
        p4d = to_ndim(Vec4f0, to_ndim(Vec3f0, p, 0f0), 1f0)
        cam_pos = view * model * p4d
        cam_pos[3] / cam_pos[4]
    end, rev=false)

    submesh = Attributes(
        model=model,
        color=color,
        shading=primitive.shading, lightposition=primitive.lightposition,
        ambient=primitive.ambient, diffuse=primitive.diffuse,
        specular=primitive.specular, shininess=primitive.shininess,
        faceculling=get(primitive, :faceculling, -10)
    )

    if !(rotations isa Vector)
        R = Makie.rotationmatrix4(to_rotation(rotations))
        submesh[:model] = model * R
    end
    scales = primitive[:markersize][]

    for i in zorder
        p = pos[i]
        if color isa AbstractVector
            submesh[:color] = color[i]
        end
        if rotations isa Vector
            R = Makie.rotationmatrix4(to_rotation(rotations[i]))
            submesh[:model] = model * R
        end
        scale = markersize isa Vector ? markersize[i] : markersize

        draw_mesh3D(
            scene, screen, submesh, mesh = m, pos = p,
            scale = scale isa Real ? Vec3f0(scale) : to_ndim(Vec3f0, scale, 1f0)
        )
    end

    return nothing
end


text:
using Dependencies
using Serialization
using Models
using HMMs
using Brobdingnag
using LinearAlgebra
using Dates
using Distributions
using Random
using Distributed

criterion( model ) = HMMs.likelihood(model.model)[end]

dir = Dependencies.defaultdir
objects = Dict{String,Any}()
for root in readdir(dir)
    file = joinpath( dir, root )
    io = open( file )
    objects[file] = deserialize( io )
    close( io )
end

statecount( ::HMMs.HMM{N} ) where {N} = N
statecount( m::Models.AbstractModel ) = statecount( Models.rootmodel( m ) )

good = statecount.( values( objects ) ) .== 3
goodmodels = collect( values( objects ) )[good];

l = length.( getfield.( Models.rootmodel.( goodmodels ), :y ) )
l .- minimum(l)

hmm = Models.rootmodel( goodmodels[3] )
C = Models.sandwich( hmm )
(eva, eve) = eigen( C )

T = real.(sqrt(C))
newhmm = deepcopy( hmm )
n = size(T,1)
p = Models.getcompressedparameters( newhmm )
Models.setcompressedparameters!( newhmm, p + T * randn(n) )

n = length(hmm.y)
sample = fill( NaN, n )
Distributions.rand!( hmm, sample )

Random.seed!(1)
N = 100
hmms = typeof(hmm)[]
for i = 1:N
    push!( hmms, rand( typeof(hmm) ) )
    HMMs.setobservations( hmms[end], sample )
end

processes = 50
addprocs( processes )
modules = [:HMMs, :Brobdingnag, :Models]
futures = Future[]
for pid in workers()
    for moduletoeval in modules
        push!( futures, remotecall( Core.eval, pid, Main, Expr(:using,Expr(:.,moduletoeval)) ) )
    end
end
for future in futures
    wait(future)
end
@time hmms = pmap( HMMs.em, hmms );

likelihoods = [HMMs.likelihood(hmm)[end] for hmm in hmms]
perm = sortperm(likelihoods)
likelihoods[perm]
hmms[perm[end]]
hmms[perm[end-1]]
hmms[perm[90]]

C = Models.sandwich( hmms[perm[end]] )
(eva, eve) = eigen(C)
B = hmms[perm[end]].basis
B * C * B'

hmm2 = deepcopy( hmm )
HMMs.setobservations( hmm2, sample )
HMMs.em( hmm2, debug=2 )
HMMs.likelihood( hmm2 )[end]
sort(likelihoods)[end-10:end]

Random.seed!(1)
N = 100
ys = [zeros(n) for i in 1:N]
[rand!( hmm, ys[i] ) for i in 1:N]
hmms = [typeof(hmm)[] for i in 1:N]
for i = 1:N
    for j = 1:N
        push!( hmms[i], rand( typeof(hmm) ) )
        HMMs.setobservations( hmms[i][end], ys[i] )
    end
end
hmms = vcat( hmms... );

@everywhere function loudfit( i, hmm )
    println( "Fitting HMM $i" )
    HMMs.em( hmm )
    println( "Done fitting HMM $i" )
    return hmm
end

fithmms = pmap( loudfit, 1:N*N, hmms );


text:
#
# Function that computes the forces and potential energy,
# using a simple quadratic energy function
#
# modifies vector f, returns u

#
# Version that receives the coordinates of the atoms only, used in equilibration
# -> does not update contaminations or count encounters
#

function uf_QUAD!(n :: Int64, x :: Array{Float64}, f :: Forces, input :: MDInput)

  resetForces!(f)

  Threads.@threads for i in 1:n-1
    it = Threads.threadid()
    for j in i+1:n
      xj = image(x[j,1] - x[i,1],input)
      yj = image(x[j,2] - x[i,2],input)
      r2 = xj^2 + yj^2
      if r2 < f.cutoff2 && r2 < f.sig2
        r = sqrt(r2)
        compute_uf_partials_QUAD!(it,f,i,j,xj,yj,r,input)
      end
    end
  end

  u = reduceForces!(f)
  return u

end

#
# Version that receives the Atoms structure, to compute the number of encounters
# -> used in the production simulation, to compute the contamination
#

function uf_QUAD!(n :: Int64, atoms :: Atoms, f :: Forces, input :: MDInput)

  x = atoms.x
  resetForces!(f)

  Threads.@threads for i in 1:n-1
    it = Threads.threadid()

    if atoms.status[i] == 2
      continue
    end

    for j in i+1:n

      if atoms.status[j] == 2
        continue
      end

      xj = image(x[j,1] - x[i,1],input)
      yj = image(x[j,2] - x[i,2],input)
      r2 = xj^2 + yj^2
      if r2 < f.cutoff2
        r = sqrt(r2)
        if r2 < f.sig2
          compute_uf_partials_QUAD!(it,f,i,j,xj,yj,r,input)
        end
        atoms.status[i], atoms.status[j], encounter = update_status(atoms.status[i],atoms.status[j],r,input)
        ipair = setipair2D(n,i,j)
        if encounter
          if ! f.encij[ipair]
            f.nenc_partial[it] = f.nenc_partial[it] + 1
          end
          f.encij[ipair] = true
        else
          f.encij[ipair] = false
        end
      end

    end

  end

  u = reduceForces!(f)
  return u, sum(f.nenc_partial)

end


text:
using Plots
using Plots.Measures
using CSV
using DataFrames
using ArgParse

pyplot()

s = ArgParseSettings()
@add_arg_table s begin
    "--scores", "-s"
        required = true
        arg_type = String
        help = "Path to the file containing the scores in tsv format (peptide protein score)."
    "--peptides", "-p"
        required = true
        arg_type = String
        help = "Path to the csv file containing the peptides data."
    "--output", "-o"
        required = true
        arg_type = String
        help = "Output file where the plots are to be saved."
end

args = parse_args(ARGS, s)

PEPTIDES_FILEPATH = args["peptides"]

SCORES_FILEPATH = args["scores"]

function load_peptide_targets(peptide_filepath::String)
	df = DataFrame(CSV.File(peptide_filepath))

	peptides = Dict{String, Array{String}}()
	for row in eachrow(df)
		peptides[row["compound_name"]] = split(row["targets"], ";")
	end

	peptides
end

function generate_curves(scores_filepath::String)
	df = DataFrame(CSV.File(scores_filepath, header=false, delim=" "))
	peptide_groups = groupby(df, :Column1)
	peptide_curves = Dict()
    println(names(peptide_groups))
	for group in peptide_groups
		peptide = group[1,:]["Column1"]
		sorted = sort(group, :Column3, rev=true)
		sorted[!, :rank] = 1:nrow(sorted)
		peptide_curves[peptide] = sorted[!,[:rank,:Column2,:Column3]]
	end
	peptide_curves
end

function plot_curve(peptide::String, curves::Dict{Any, Any}, targets::Dict{String, Array{String}})
	df = curves[peptide]
	non_targets = filter(x -> x["Column2"] ∉ targets[peptide], df)
	targets = filter(x -> x["Column2"] ∈ targets[peptide], df)
	scatter(
		non_targets[!, :rank],
		non_targets[!, :Column3],
	  	title="O2A Curve ($peptide)",
		label="Non-target(s)",
		xlabel="Rank", ylabel="Score",
        margin=10mm
       )
	scatter!(
		targets[!, :rank],
		targets[!, :Column3],
		markersize=6,
		color="red", label="Target(s)")
end

println("Loading the peptide data...")
targets = load_peptide_targets(PEPTIDES_FILEPATH)

println("Generating the O2A curves...")
curves = generate_curves(SCORES_FILEPATH)

println("Plotting the curves...")
plots = [plot_curve(peptide, curves, targets) for peptide in keys(targets)]
plot(plots..., layout = length(plots), legend = false, size = (2000,2000), margin=20mm)

println("Saving the plots to $(args["output"])...")
savefig(args["output"])

max_stars_repo_name: GreenCUBIC/PeptideScreening

max_issues_repo_name: GreenCUBIC/PeptideScreening

max_forks_repo_name: GreenCUBIC/PeptideScreening


text:
using JSON3
using MacroTools
using JuliaFormatter: format_text, format_file

import SPIRV_Headers_jll

const include_dir = joinpath(SPIRV_Headers_jll.artifact_dir, "include", "spirv", "unified1")

include("formatting.jl")
include("grammar.jl")

mkpath(src_dir("generated"))

generate() && format()


text:
module wembed
using HDF5
using Base.Collections

# core.jl
  export init, bin2hdf5, getsimilarity, getvector, wordnearest, odd1out, nsimilarity

  include("core.jl")
end # module


text:
struct EPOCHSimulation{P,B,C,PC}
    dir::String
    files::Vector{SDFFile{P,B,C,PC}}
    param::P
    cache::C
    particle_cache::PC
end

function read_simulation(dir; field_cache_size=6, particle_cache_size=2, kwargs...)
    paths = filter(f->endswith(f, ".sdf"), readdir(dir))
    !issorted(paths) && sort!(paths)
    if isempty(paths)
        error("no .sdf files found in directory $dir")
        return nothing
    end

    input_deck = joinpath(dir, "input.deck")
    if isfile(input_deck)
        p = parse_input(input_deck)
    else
        @warn "No input.deck found in $dir."
        p = missing
    end

    # We assume that the particle_cache has the same dimensionality throughout
    # all the simulation, so we ca use the first file to get the dimensionality
    # and element type

    N, T = get_data_description(joinpath(dir, first(paths)))
    @debug "First file gave N = $N and T = $T"

    field_cache = LRU{Tuple{String,AbstractBlockHeader},AbstractArray}(maxsize=field_cache_size, kwargs...)
    particle_cache =  LRU{Tuple{String,AbstractBlockHeader},NTuple{N,Vector{T}}}(maxsize=particle_cache_size)

    files = read_file.(joinpath.((dir,), paths), (Ref(p),), (Ref(field_cache),), (Ref(particle_cache),))

    EPOCHSimulation(dir, files, p, field_cache, particle_cache)
end

# Indexing
Base.getindex(sim::EPOCHSimulation, i::Int) = sim.files[i]
Base.firstindex(sim::EPOCHSimulation) = firstindex(sim.files)
Base.lastindex(sim::EPOCHSimulation) = lastindex(sim.files)

# Iteration
Base.iterate(sim::EPOCHSimulation, state...) = iterate(sim.files, state...)
Base.eltype(::Type{EPOCHSimulation}) = SDFFile
Base.length(sim::EPOCHSimulation) = length(sim.files)
Base.size(sim::EPOCHSimulation, args...) = size(sim.files, args...)

# Utils

Base.ndims(sim::EPOCHSimulation) = ndims(first(sim))
Base.haskey(sim::EPOCHSimulation, key) = haskey(sim.param, key)
Base.haskey(sim::EPOCHSimulation, block, key) = haskey(sim.param[block], key)

# show
function Base.show(io::IO, ::MIME"text/plain", sim::EPOCHSimulation)
    first_file = first(sim)
    code_name = first_file.header.code_name
    n = length(sim)
    t₀ = si_round(get_time(first_file))
    t = si_round(get_time(last(sim)))
    fc = Base.summarysize(sim.cache)
    pc = Base.summarysize(sim.particle_cache)
    c = Base.format_bytes(fc + pc)
    description = "$code_name simulation with $n files from " * t₀ * " to " * t * ".\n" *
    "Cache size: " * c
    print(io, description)
end


text:
module TestRecyclers

using Random: rand!
using Recyclers
using Test

using ..Utils

function check_seq_untyped(R, factory)
    T = typeof(factory())
    recycler = R(factory)
    @test recycler isa R{T}
    check_seq(T, recycler)
end

function check_seq_typed(R, factory)
    T = typeof(factory())
    check_seq(T, R{T}(factory))
end

function check_seq(T, recycler)
    @test eltype(recycler) === T
    @test Recyclers.get!(recycler) isa T
    @test Recyclers.get!(recycler) isa T
    object = Recyclers.get!(recycler)
    @test Recyclers.recycle!(recycler, object)
    @test Recyclers.get!(recycler) === object
end

zeros3() = zeros(Int, 3)

function test_seq_untyped()
    @testset "$(nameof(R))" for R in Utils.recycler_types()
        check_seq_untyped(R, zeros3)
    end
end

function test_seq_typed()
    @testset "$(nameof(R))" for R in Utils.recycler_types()
        check_seq_typed(R, zeros3)
    end
end

function check_concurrency(recycler; ntrials = 2^10, ntasks = Threads.nthreads())
    failures = [Ref(0) for _ in 1:ntasks]
    @sync for fs in failures
        Threads.@spawn begin
            tmp = Recyclers.get!(recycler)
            for _ in 1:ntrials
                Recyclers.recycling!(recycler) do xs
                    rand!(tmp)
                    copyto!(xs, tmp)
                    rand(Bool) && yield()
                    fs[] += !isequal(xs, tmp)
                end
            end
        end
    end
    @test sum(getindex, failures) == 0
end

# Since `check_concurrency` can be time-confusing, separate test cases to functions so that
# `TestFunctionRunner` can terminate the test at a more granular unit.
for R in Utils.recycler_types(), oversubscribe in [1, 3]
    @eval function $(Symbol(:test_concurrency_, nameof(R), :_, Symbol(oversubscribe)))(
        ntrials = 2^10,
    )
        check_concurrency($R(zeros3); ntasks = Threads.nthreads() * $oversubscribe)
    end
end

function test_reinit_sharded()
    @testset "$(nameof(R))" for R in Utils.recycler_types()
        recycler = R(zeros3)
        @test Recyclers.unsafe_destroy!(recycler) === recycler
        @test Recyclers.unsafe_init!(recycler) === recycler
        if recycler isa Recyclers.ShardedRecycler
            @test length(recycler.cache.pools) == Threads.nthreads()
        end
    end
end

end  # module


text: precompile(Tuple{typeof(Base.peek), Base.IOStream})


text:
# Test LZ78 Predictor

@testset "LZ78" begin

p       = LZ78{Char}();
data    = [ 'a', 'a', 'a', 'b', 'a',
            'b', 'b', 'b', 'b', 'b',
            'a', 'a', 'b', 'c', 'c',
            'd', 'd', 'c', 'b', 'a',
            'a', 'a', 'a' ];

# # Create model
@testset "Online Prediction" begin
    # Create model
    for symbol ∈ data
        add!( p, symbol );
        predictions     = predict( p );
        @test mapreduce(sym->predictions[sym],+,keys(predictions)) <= 1.00001;     # :/ not adding upto 1.0
    end
end

# Test model
test_model = Dict{Vector{Char},Int}(
    ['a'] => 5,
    ['a','a'] => 2,
    ['a','a','a'] => 1,
    ['a','b'] => 2,
    ['a','b','c'] => 1,
    ['b'] => 4,
    ['b','a'] => 1,
    ['b','b'] => 2,
    ['b','b','a'] => 1,
    ['c'] => 1,
    ['d'] => 2,
    ['d','c'] => 1
);

@testset "Model Tests" begin
    # Verify Model
    @testset "Values" begin
        @test p.model.value == 23;
        for node ∈ keys( test_model )
            @test p.model[node] == test_model[node];
        end
    end

    # Verify no other nodes exists
    @testset "Nodes" begin
        nodes   = keys( test_model )
        for node in keys( p.model )
            @test node ∈ nodes;
        end
    end

    @testset "Symbol Prediction" begin
        sym_probability = predict( p );
        @test sym_probability['a']  ≈  89/120;
        @test sym_probability['b']  ≈  28/120;
        @test sym_probability['c']  ≈   1/120;
        @test sym_probability['d']  ≈   2/120;
    end

    @testset "Model Size" begin
        @test size( p ) == 13;
    end
end

@testset "Best Symbol" begin
    # Reusing p
    @test get_best_symbol( p ) == 'a'
    # For untrained model
    @test get_best_symbol( LZ78{Int}() ) == nothing
    @test get_best_symbol( LZ78{Char}() ) == nothing
    @test get_best_symbol( LZ78{String}() ) == nothing
end

@testset "Info tests" begin
    @test info_string( p ) == "LZ78";
    @test unique_string( p ) == "LZ78";
end

@testset "Random Sequence test" begin
    for idx = 1:10
        learnt_model    = LZ78{Int}();
        @test @test_nothrow for i = 1:1000
            symbol = trunc(Int,10*rand());
            add!( learnt_model, symbol );
            predict( learnt_model );
        end
    end
end

end

nothing


text:
#    Copyright (c) 2021 J.A. Duffek
#    Copyright (c) 2000 D.M. Spink
#
#    This program is free software; you can redistribute it and/or modify
#    it under the terms of the GNU General Public License as published by
#    the Free Software Foundation; either version 2 of the License, or
#    (at your option) any later version.
#
#    This program is distributed in the hope that it will be useful,
#    but WITHOUT ANY WARRANTY; without even the implied warranty of
#    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#    GNU General Public License for more details.
#
#    You should have received a copy of the GNU General Public License
#    along with this program; if not, see <http://www.gnu.org/licenses/>.

"""
    vecangle(num::Matrix{F},den::Matrix{F}
            )::Matrix{F} where {F<:AbstractFloat}

or

    vecangle(num::F,den::F)::Vector{F} where {F<:AbstractFloat}

An alternative to atan, returning an arctangent in the range `0` to `2*π`.

The components of the vector ang are the arctangent of the corresponding enties
of `num./den`. This function is an alternative for atan, returning an angle in
the range `0` to `2*π`.

# Arguments:
- `num`: Numerator, vector of size `(1,nv)`.
- `den`: Denominator, vector of size `(1,nv)`.

# Output:
- `ang`: Arctangents, row vector of angles.

# Examples:
```julia
julia> ang = vecmag2(num,dum)
```

Find the `atan(1.2,2.0)` and `atan(1.5,3.4)` using `vecangle`.
```julia
julia> ang = vecangle([1.2 1.5], [2.0 3.4])
```
"""
function vecangle(num::Matrix{F},den::Matrix{F}
                 )::Matrix{F} where {F<:AbstractFloat}
ang = atan.(num,den);
ang[ang .< 0.0] .+= 2*π;
return ang;
end # vecangle

# overloaded function for float input to get a float instead of a vector
function vecangle(num::F,den::F)::F where {F<:AbstractFloat}
ang = atan(num,den);
if ang < 0.0
  ang += 2*π;
end # if
return ang;
end # vecangle


text:
"""
The ER sampler consists of:
    - N, the set of nodes
    - k, Integer, hyperedge size
In the ER model, all edges have the same probability of placement.
You may use the function ER_sampler below to generate the sampler.
"""
struct ER_sampler{N, K} <: Random.Sampler{NTuple{K, eltype(N)}}
    nodes::N
end

"""
Generate the ER_sampler
2 parameters:
- n :: Integer, node count
- k :: Integer, hyperedge size
"""
function ER_sampler(n::Integer, k::Integer)
    ER_sampler{Base.OneTo{typeof(n)}, k}(Base.OneTo(n))
end

function Base.rand(rng::AbstractRNG, s::ER_sampler{T, K}) where {T, K}
    ntuple(i -> rand(rng, s.nodes), Val{K}())
end

max_stars_repo_name: LilithHafner/FBD

max_issues_repo_name: LilithHafner/FBD

max_forks_repo_name: LilithHafner/FBD


text:
using StackOverflow

"""
    searcherror(intitle::String)

Adding the ability to search for an error through StackOverflow.jl while in the Julia Terminal.
"""
function searcherror(intitle::String)
    data = replace(intitle, " " => "%20")
    r = HTTP.request("GET", "https://api.stackexchange.com/2.2/search?order=desc&sort=activity&intitle=$(data)&site=stackoverflow")
    json_obj = StackOverflow.convert_HTTP_Response_To_JSON(r)
    makequestionsarray(json_obj)
end


text:
"""
    fixPoint!(model::CPModel, new_constraints=nothing)

Run the fix-point algorithm. Will prune the domain of every variable of `model` as much as possible, using its constraints.
Return a tuple with a boolean corresponding to the feasibility and a `CPModification` object, containing all the pruned domains.

# Arguments
- `model::CPModel`: the model you want to apply the algorithm on.
- `new_constraints::Union{Array{Constraint}, Nothing}`: if this is filled with a set of constraints,
only those will be propagated in the first place.
"""
function fixPoint!(model::CPModel, new_constraints::Union{Array{Constraint}, Nothing}=nothing, prunedDomains::Union{CPModification,Nothing}=nothing)
    toPropagate = Set{Constraint}()

    if isnothing(prunedDomains)
        prunedDomains = CPModification()
    end

    # If we did not specify the second argument, it is the beginning so we propagate every constraint
    if isnothing(new_constraints)
        addToPropagate!(toPropagate, model.constraints)
    else
        addToPropagate!(toPropagate, new_constraints)
    end

    # Dealing with the objective
    if !isnothing(model.objectiveBound)
        prunedObj = removeAbove!(model.objective.domain, model.objectiveBound)
        if isempty(model.objective.domain)
            return false, prunedDomains
        end

        if !isempty(prunedObj)
            addToPrunedDomains!(prunedDomains, model.objective, prunedObj)
            triggerDomainChange!(toPropagate, model.objective)
        end
    end

    while !isempty(toPropagate)
        constraint = pop!(toPropagate)
        if !propagate!(constraint, toPropagate, prunedDomains)
            triggerInfeasible!(constraint, model)
            return false, prunedDomains
        end
    end

    return true, prunedDomains
end


text:
# Polynomial Normal Form

"""
    labels!(dict, t)

Find all terms that are not + and * and replace them
with a symbol, store the symbol => term mapping in `dict`.
"""
function labels! end

# Turn a Term into a multivariate polynomial
function labels!(dicts, t::Sym,  variable_type::Type)
    sym2term, term2sym = dicts
    if !haskey(term2sym, t)
        sym2term[t] = t
        term2sym[t] = t
    end
    return t
end

function labels!(dicts, t, variable_type::Type)
    if t isa Number
        return t
    elseif istree(t) && (operation(t) == (*) || operation(t) == (+) || operation(t) == (-))
        tt = arguments(t)
        return similarterm(t, operation(t), map(x->labels!(dicts, x, variable_type), tt), symtype(t))
    elseif istree(t) && operation(t) == (^) && length(arguments(t)) > 1 && isnonnegint(arguments(t)[2])
        return similarterm(t, operation(t), map(x->labels!(dicts, x, variable_type), arguments(t)), symtype(t))
    else
        sym2term, term2sym = dicts
        if haskey(term2sym, t)
            return term2sym[t]
        end
        if istree(t)
            tt = arguments(t)
            sym = Sym{symtype(t)}(gensym(nameof(operation(t))))
            dicts2 = _dicts(dicts[2])
            sym2term[sym] = similarterm(t, operation(t),
                                        map(x->to_mpoly(x, variable_type, dicts)[1], arguments(t)),
                                        symtype(t))
        else
            sym = Sym{symtype(t)}(gensym("literal"))
            sym2term[sym] = t
        end

        term2sym[t] = sym

        return sym
    end
end

ismpoly(x) = x isa MP.AbstractPolynomialLike || x isa Number
isnonnegint(x) = x isa Integer && x >= 0

_dicts(t2s=OrderedDict{Any, Sym}()) = (OrderedDict{Sym, Any}(), t2s)

let
    mpoly_preprocess = [@rule(identity(~x) => ~x)
                        @rule(zero(~x) => 0)
                        @rule(one(~x) => 1)]

    simterm(x, f, args;metadata=nothing) = similarterm(x,f,args, symtype(x); metadata=metadata)
    mpoly_rules = [@rule(~x::ismpoly - ~y::ismpoly => ~x + -1 * (~y))
                   @rule(-(~x) => -1 * ~x)
                   @acrule(~x::ismpoly + ~y::ismpoly => ~x + ~y)
                   @rule(+(~x) => ~x)
                   @acrule(~x::ismpoly * ~y::ismpoly => ~x * ~y)
                   @rule(*(~x) => ~x)
                   @rule((~x::ismpoly)^(~a::isnonnegint) => (~x)^(~a))]
    global const MPOLY_CLEANUP = Fixpoint(Postwalk(PassThrough(RestartedChain(mpoly_preprocess)), similarterm=simterm))
    MPOLY_MAKER = Fixpoint(Postwalk(PassThrough(RestartedChain(mpoly_rules)), similarterm=simterm))

    global to_mpoly
    function to_mpoly(t, variable_type::Type=DynamicPolynomials.PolyVar{true}, dicts=_dicts())
        # term2sym is only used to assign the same
        # symbol for the same term -- in other words,
        # it does common subexpression elimination
        t = MPOLY_CLEANUP(t)
        sym2term, term2sym = dicts
        labeled = labels!((sym2term, term2sym), t, variable_type)

        if isempty(sym2term)
            return MPOLY_MAKER(labeled), Dict{Sym,Any}()
        end

        ks = sort(collect(keys(sym2term)), lt=<ₑ)
        vars = MP.similarvariable.(variable_type, nameof.(ks))

        replace_with_poly = Dict{Sym,eltype(vars)}(zip(ks, vars))
        t_poly = substitute(labeled, replace_with_poly, fold=false)
        MPOLY_MAKER(t_poly), sym2term
    end
end

function to_term(reference, x, dict)
    syms = Dict(zip(string.(nameof.(keys(dict))), keys(dict)))
    dict = copy(dict)
    for (k, v) in dict
        dict[k] = _to_term(reference, v, dict, syms)
    end
    return _to_term(reference, x, dict, syms)
    #return substitute(t, dict, fold=false)
end

_to_term(reference, x::Number, dict, syms) = x
_to_term(reference, var::MP.AbstractVariable, dict, syms) = substitute(syms[MP.name(var)], dict, fold=false)
function _to_term(reference, mono::MP.AbstractMonomialLike, dict, syms)
    monics = [
        begin
            t = _to_term(reference, var, dict, syms)
            exp == 1 ? t : t^exp
        end
        for (var, exp) in MP.powers(mono) if !iszero(exp)
    ]
    if length(monics) == 1
        return monics[1]
    elseif isempty(monics)
        return 1
    else
        return similarterm(reference, *, monics, symtype(reference))
    end
end

function _to_term(reference, term::MP.AbstractTermLike, dict, syms)
    coef = MP.coefficient(term)
    mono = _to_term(reference, MP.monomial(term), dict, syms)
    if isone(coef)
        return mono
    else
        return MP.coefficient(term) * mono
    end
end

function _to_term(reference, x::MP.AbstractPolynomialLike, dict, syms)
    if MP.nterms(x) == 0
        return 0
    elseif MP.nterms(x) == 1
        return _to_term(reference, first(MP.terms(x)), dict, syms)
    else
        terms = map(MP.terms(x)) do term
            _to_term(reference, term, dict, syms)
        end
        return similarterm(reference, +, terms, symtype(reference))
    end
end

function _to_term(reference, x, dict, vars)
    if istree(x)
        t = similarterm(x, operation(x), _to_term.((reference,), arguments(x), (dict,), (vars,)), symtype(x))
    else
        if haskey(dict, x)
            return dict[x]
        else
            return x
        end
    end
end

<ₑ(a::MP.AbstractPolynomialLike, b::MP.AbstractPolynomialLike) = false

"""
    expand(expr, variable_type::Type=DynamicPolynomials.PolyVar{true})

Expand expressions by distributing multiplication over addition, e.g.,
`a*(b+c)` becomes `ab+ac`.

`expand` uses replace symbols and non-algebraic expressions by variables of type
`variable_type` to compute the distribution using a specialized sparse
multivariate polynomials implementation.
`variable_type` can be any subtype of `MultivariatePolynomials.AbstractVariable`.
"""
function expand(expr, variable_type::Type=DynamicPolynomials.PolyVar{true})
    to_term(expr, to_mpoly(expr, variable_type)...)
end

## Hack to fix https://github.com/JuliaAlgebra/MultivariatePolynomials.jl/issues/169

Base.promote_rule(::Type{S}, ::Type{T}) where {S<:Symbolic, T<:MP.AbstractPolynomialLike}= Any
Base.promote_rule(::Type{T}, ::Type{S}) where {S<:Symbolic, T<:MP.AbstractPolynomialLike}= Any

Base.@deprecate polynormalize(x) expand(x)


text:
module EchoREPL

import REPL
import REPL.LineEdit

do_cmd(repl, input) = println(REPL.outstream(repl), input)

function create_mode(repl, main_mode)
    echo_mode = LineEdit.Prompt("echo> ";
        prompt_prefix = Base.text_colors[:cyan],
        prompt_suffix = "",
        # complete = nothing,
        sticky = false)
    echo_mode.repl = repl
    hp = main_mode.hist
    hp.mode_mapping[:echo] = echo_mode
    echo_mode.hist = hp

    search_prompt, skeymap = LineEdit.setup_search_keymap(hp)
    prefix_prompt, prefix_keymap = LineEdit.setup_prefix_keymap(hp, echo_mode)

    echo_mode.on_done = (s, buf, ok) -> begin
        ok || return REPL.transition(s, :abort)
        input = String(take!(buf))
        REPL.reset(repl)
        do_cmd(repl, input)
        REPL.prepare_next(repl)
        REPL.reset_state(s)
        s.current_mode.sticky || REPL.transition(s, main_mode)
    end

    mk = REPL.mode_keymap(main_mode)

    # shell_mode = nothing
    # for mode in Base.active_repl.interface.modes
    #     if mode isa LineEdit.Prompt
    #         mode.prompt == "shell> " && (shell_mode = mode)
    #     end
    # end

    # repl_keymap = Dict()
    # if shell_mode != nothing
    #     repl_keymap[';'] = function (s, o...)
    #         if isempty(s) || position(LineEdit.buffer(s)) == 0
    #             buf = copy(LineEdit.buffer(s))
    #             LineEdit.transition(s, shell_mode) do
    #                 LineEdit.state(s, shell_mode).input_buffer = buf
    #             end
    #         else
    #             LineEdit.edit_insert(s, ';')
    #         end
    #     end
    # end

    b = Dict{Any, Any}[
        skeymap, #=repl_keymap,=# mk, prefix_keymap, LineEdit.history_keymap,
        LineEdit.default_keymap, LineEdit.escape_defaults
    ]
    echo_mode.keymap_dict = LineEdit.keymap(b)
    return echo_mode
end

function repl_init(repl)
    main_mode = repl.interface.modes[1]
    echo_mode = create_mode(repl, main_mode)
    push!(repl.interface.modes, echo_mode)
    org_action = main_mode.keymap_dict['\x01']
    keymap = Dict{Any, Any}(
        # "^A" => function (s, args...)
        '\x01' => function (s, args...)
            # println("^A")
            if isempty(s) || position(LineEdit.buffer(s)) == 0
                buf = copy(LineEdit.buffer(s))
                LineEdit.transition(s, echo_mode) do
                    LineEdit.state(s, echo_mode).input_buffer = buf
                end
            else
                # LineEdit.edit_insert(s, ']')
                org_action(s, args...)
            end
        end
    )
    main_mode.keymap_dict = LineEdit.keymap_merge(main_mode.keymap_dict, keymap)
    return
end

function __init__()
    if isdefined(Base, :active_repl)
        repl_init(Base.active_repl)
    else
        atreplinit() do repl
            if isinteractive() && repl isa REPL.LineEditREPL
                isdefined(repl, :interface) || (repl.interface = REPL.setup_interface(repl))
                repl_init(repl)
            end
        end
    end
end

end # module


text:
using ObjectDetector, FileIO

yolomod = YOLO.v3_608_COCO(batch=1, silent=true) # Load the YOLOv3-tiny model pretrained on COCO, with a batch size of 1

batch = emptybatch(yolomod) # Create a batch object. Automatically uses the GPU if available

img = load("images.jpeg")

batch[:,:,:,1], padding = prepareImage(img, yolomod) # Send resized image to the batch

res = yolomod(batch, detectThresh=0.5, overlapThresh=0.8) # Run the model on the length-1 batch

imgBoxes = drawBoxes(img, yolomod, padding, res)
save("result.png", imgBoxes)

max_stars_repo_name: ashwani-rathee/Sigil

max_issues_repo_name: ashwani-rathee/Sigil

max_forks_repo_name: ashwani-rathee/Sigil


text:
# @everywhere f(s,count)=(println("process id = $(myid()) s = $s count = $count");repeat(s,count))
# pmap((a1,a2)->f(a1,a2),{"a","b","c"},{2,1,3})

@everywhere type Foo
	r
end

@everywhere function update!(el,n)
	el.r = -1 * n
end

x = [Foo(2*i) for i=1:1000]

d = [i=>i for i=1:1000]
samples_pmap = pmap(update!,x,d)
samples = vcat(samples_pmap[1],samples_pmap[2],samples_pmap[3],samples_pmap[4])
print(x)


text:
# This file includes testset for DAESystem

@testset "DAESystemTestSet" begin
    @info "Running DAESystemTestSet ..."

    # DAESystem construction
    function sfunc(out, dx, x, u, t)
        out[1] = -0.04 * x[1] + 1e4 * x[2] * x[3] - dx[1]
        out[2] = 0.04 * x[1] - 3e7 * x[2]^2 - 1e4 * x[2] * x[3] - dx[2]
        out[3] = x[1] + x[2] + x[3] - 1.0
    end
    ofunc(x, u, t) = x
    state = [1., 0., 0.]
    stateder = [-0.04, 0.04, 0.]
    differential_vars = [true, true, false]
    ds = DAESystem(righthandside=sfunc, readout=ofunc, state=state, input=nothing, output=Outport(3),
        stateder=stateder, diffvars=differential_vars)
    @test typeof(ds.trigger) == Inpin{Float64}
    @test typeof(ds.handshake) == Outpin{Bool}
    @test ds.input === nothing
    @test typeof(ds.output) == Outport{Outpin{Float64}}
    @test length(ds.output) == 3
    @test ds.state == state
    @test ds.t == 0.
    @test ds.integrator.sol.prob.p === nothing

    # Driving DAESystem
    iport = Inport(3)
    trg = Outpin()
    hnd = Inpin{Bool}()
    connect!(ds.output, iport)
    connect!(trg, ds.trigger)
    connect!(ds.handshake, hnd)
    tsk = launch(ds)
    tsk2 = @async while true
        all(take!(iport) .=== NaN) && break
    end
    for t in 1 : 10
        put!(trg, t)
        take!(hnd)
        @test ds.t == t
        @test [read(pin.link.buffer) for pin in iport] == ds.state
    end
    put!(trg, NaN)
    sleep(0.1)
    @test istaskdone(tsk)
    put!(ds.output, NaN * ones(length(ds.output)))
    sleep(0.1)
    @test istaskdone(tsk2)

    # DAESystem with inputs
    function sfunc2(out, dx, x, u, t)
        out[1] = -0.04 * x[1] + 1e4 * x[2] * x[3] - dx[1] + u[1](t)
        out[2] = 0.04 * x[1] - 3e7 * x[2]^2 - 1e4 * x[2] * x[3] - dx[2] + u[2](t)
        out[3] = x[1] + x[2] + x[3] - 1.0
    end
    ofunc2(x, u, t) = x
    state = [1., 0., 0.]
    stateder = [-0.04, 0.04, 0.]
    differential_vars = [true, true, false]
    ds = DAESystem(righthandside=sfunc2, readout=ofunc2, state=state, input=Inport(2), output=Outport(3),
        stateder=stateder, diffvars=differential_vars)
    @test typeof(ds.trigger) == Inpin{Float64}
    @test typeof(ds.handshake) == Outpin{Bool}
    @test typeof(ds.input) <: Inport
    @test typeof(ds.output) == Outport{Outpin{Float64}}
    @test length(ds.input) == 2
    @test length(ds.output) == 3
    @test ds.state == state
    @test ds.t == 0.
    @test typeof(ds.integrator.sol.prob.p) <: Interpolant
    @test size(ds.integrator.sol.prob.p.timebuf) == (3,)
    @test size(ds.integrator.sol.prob.p.databuf) == (2,3)

    # Driving DAESystem with input
    oport = Outport(2)
    iport = Inport(3)
    trg = Outpin()
    hnd = Inpin{Bool}()
    connect!(oport, ds.input)
    connect!(ds.output, iport)
    connect!(trg, ds.trigger)
    connect!(ds.handshake, hnd)
    tsk = launch(ds)
    tsk2 = @async while true
        all(take!(iport) .=== NaN) && break
    end
    for t in 1 : 10
        put!(trg, t)
        put!(oport, ones(2) * t)
        take!(hnd)
        @test ds.t == t
        @test [read(pin.link.buffer) for pin in iport] == ds.state
    end
    put!(trg, NaN)
    sleep(0.1)
    @test istaskdone(tsk)
    put!(ds.output, NaN * ones(length(ds.output)))
    sleep(0.1)
    @test istaskdone(tsk2)

    # Test defining new DAESystems
    # Type mest be mutable
    @test_throws Exception @eval @def_dae_system struct DAESystem{RH, RO, ST, IP, OP}
        righthandside::RH
        readout::RO
        state::ST
        stateder::ST
        diffvars::Vector{Bool}
        input::IP
        output::OP
    end

    # The type must be a subtype of AbstractDAESystem
    @test_throws Exception @eval @def_dae_system mutable struct DAESystem{RH, RO, ST, IP, OP}
        righthandside::RH
        readout::RO
        state::ST
        stateder::ST
        diffvars::Vector{Bool}
        input::IP
        output::OP
    end

    # The type must be a subtype of AbstractDAESystem
    @test_throws Exception @eval @def_dae_system mutable struct DAESystem{RH, RO, ST, IP, OP} <: MyDummyAbstractDAESystem
        righthandside::RH
        readout::RO
        state::ST
        stateder::ST
        diffvars::Vector{Bool}
        input::IP
        output::OP
    end

    @info "Done DAESystemTestSet."
end # testset


text:
@doc raw"""
    SubGradientProblem <: Problem

A structure to store information about a subgradient based optimization problem

# Fields
* `manifold` – a [Manifold](https://juliamanifolds.github.io/Manifolds.jl/stable/interface.html#ManifoldsBase.Manifold)
* `cost` – the function $F$ to be minimized
* `subgradient` – a function returning a subgradient $\partial F$ of $F$

# Constructor

    SubGradientProblem(M, f, ∂f)

Generate the [`Problem`] for a subgradient problem, i.e. a function `f` on the
manifold `M` and a function `∂f` that returns an element from the subdifferential
at a point.
"""
struct SubGradientProblem{mT<:Manifold,TCost,TSubgradient} <: Problem
    M::mT
    cost::TCost
    subgradient::TSubgradient
end
"""
    get_subgradient(p,x)

Evaluate the (sub)gradient of a [`SubGradientProblem`](@ref)` p` at the point `x`.
"""
function get_subgradient(p::SubGradientProblem, x)
    return p.subgradient(x)
end
"""
    SubGradientMethodOptions <: Options
stories option values for a [`subgradient_method`](@ref) solver

# Fields
* `retraction_method` – the retration to use within
* `stepsize` – a [`Stepsize`](@ref)
* `stop` – a [`StoppingCriterion`](@ref)
* `x` – (initial or current) value the algorithm is at
* `x_optimal` – optimal value
* `∂` the current element from the possivle subgradients at `x` that is used
"""
mutable struct SubGradientMethodOptions{TRetract<:AbstractRetractionMethod,TStepsize,P,T} <:
               Options where {P,T}
    retraction_method::TRetract
    stepsize::TStepsize
    stop::StoppingCriterion
    x::P
    x_optimal::P
    ∂::T
    function SubGradientMethodOptions(
        M::TM,
        x::P,
        sC::StoppingCriterion,
        s::Stepsize,
        retraction_method=ExponentialRetraction(),
    ) where {TM<:Manifold,P}
        return new{typeof(retraction_method),typeof(s),P,typeof(zero_tangent_vector(M, x))}(
            retraction_method, s, sC, x, x
        )
    end
end


text:
#This Script evaluates multiple Algorithms for the desired Environment for training comparison!
include("src/init.jl")

"""
REINFORCE
"""
function reinforce(params::Params, θa::Any, θv::Any, res::Array{Float64,1}, prec::Array{Float64,1}, all_traj::Traj_Data, best_net::Any)
    #Init Part:
    best_net = copy(θa)
    min_cost = Inf
    #------------------------
    #Initialize the Policy
    ĝ_a = copy(θa)
    opt_para_a = optimizers(θa,   Adam; lr = params.lr_pl)

    #Initialize the Value Function
    ĝ_v = copy(θv)
    opt_para_v = optimizers(θv,   Adam; lr = params.lr_vf)

    #Init Progress bar
    progress = Progress(params.n_iter,1)

    # Init rms
    rms_all = init_rms(params.s_init)
    rms_loc = init_rms(params.s_init)

    a_t_old = all_traj.traj[1].a_t
    weights_old = zeros(length(net_to_vec(θa)),2)

    for i_iter = 1:params.n_iter
        ĝ_a *= 0
        ĝ_v *= 0

        #Save data in intervalls
        if floor(i_iter/50) == ceil(i_iter/50)
            #make_video(all_traj,params, i_iter)
            #save_data(res, prec, θa, θv, all_traj, i_iter)
        end

        #-----------------------------------------------------------------------
        #Solve it in Parallel!
        res[i_iter], prec[i_iter] = parallel_episodes(params,θa,θv,all_traj,rms_all)

        #Apply simple baseline
        simple_baseline!(all_traj::Traj_Data, params::Params)

        #Reorganice Data in minibatches
        mini_batch = build_minibatches(all_traj,params);

        #Apply the gradient
        ĝ_a = parallel_gradient(params, θa, mini_batch)
        Knet.update!(θa, ĝ_a, opt_para_a)

        if min_cost > res[i_iter]
            min_cost = res[i_iter]
            best_net = copy(θa)
            print(" new_best: ", round(min_cost, digits = 3))
        end

        #Visualization_stuff
        #-----------------------------------------------------------------------
        ProgressMeter.update!(progress,i_iter) #update progress bar
        if floor(i_iter/50) == ceil(i_iter/50)
            a_t_old ,weights_old = plot_stuff(params, all_traj, a_t_old,
            weights_old, i_iter,θa, θv, res, prec/params.n_traj, ĝ_a, ĝ_v)
        end
        #-----------------------------------------------------------------------

        #clean the workspace from junk
        @everywhere GC.gc(true)
    end
    return res
end

"""
ADVANTAGE ACTOR CRITIC
"""
function advantage_actor_critic(params::Params, θa::Any, θv::Any, res::Array{Float64,1}, prec::Array{Float64,1}, all_traj::Traj_Data, best_net::Any)
    #Init Part:
    best_net = copy(θa)
    min_cost = Inf
    #------------------------
    #Initialize the Policy
    ĝ_a = copy(θa)
    opt_para_a = optimizers(θa,   Adam; lr = params.lr_pl)

    #Initialize the Value Function
    ĝ_v = copy(θv)
    opt_para_v = optimizers(θv,   Adam; lr = params.lr_vf)

    #Init Progress bar
    progress = Progress(params.n_iter,1)

    # Init rms
    rms_all = init_rms(params.s_init)
    rms_loc = init_rms(params.s_init)

    kl_means = zeros(params.n_iter)

    #saving old data to see differences
    a_t_old = all_traj.traj[1].a_t
    weights_old = zeros(length(net_to_vec(θa)),2)

    for i_iter = 1:params.n_iter
        ĝ_a *= 0
        ĝ_v *= 0

        #execute the parallel_episodes
        res[i_iter], prec[i_iter] = parallel_episodes(params,θa,θv,all_traj,rms_all)

        #Reorganice Data in minibatches
        mini_batch = build_minibatches(all_traj,params)

        #First refit the value function for new data batch:
        θv = refit_VF(params, θv, mini_batch, opt_para_v)

        #Now apply the gae for the new data based on updated value function
        gae_all!(all_traj,θv,params)

        #Reorganice Data in minibatches as advantages have changed...
        mini_batch = build_minibatches(all_traj,params)

        #Get Policy Update
        ĝ_a = parallel_gradient(params, θa, mini_batch)
        Knet.update!(θa, ĝ_a, opt_para_a)

        #Get KL-div of last update step
        kl_means[i_iter] = kl_all(θa,all_traj::Traj_Data,params::Params)

        if min_cost > res[i_iter]
            min_cost = res[i_iter]
            best_net = copy(θa)
            print(" new_best: ", round(min_cost, digits = 3))
        end

        #Visualization_stuff
        #-----------------------------------------------------------------------
        ProgressMeter.update!(progress,i_iter) #update progress bar
        if floor(i_iter/50) == ceil(i_iter/50)
            a_t_old ,weights_old = plot_stuff(params, all_traj, a_t_old,
            weights_old, i_iter,θa, θv, res, prec/params.n_traj, ĝ_a, ĝ_v)
        end
        #-----------------------------------------------------------------------

        @everywhere GC.gc(true)
    end
    return res, kl_means
end

"""
TNPG
"""
function tnpg(params::Params, θa::Any, θv::Any, res::Array{Float64,1}, prec::Array{Float64,1}, all_traj::Traj_Data, best_net::Any)
    #Init Part:
    best_net = copy(θa)
    min_cost = Inf
    #------------------------
    #Initialize the Policy
    ĝ_a = copy(θa)
    opt_para_a = optimizers(θa,   Adam; lr = params.lr_pl)

    #Initialize the Value Function
    ĝ_v = copy(θv)
    opt_para_v = optimizers(θv,   Adam; lr = params.lr_vf)

    #Init Progress bar
    progress = Progress(params.n_iter,1)

    # Init rms
    rms_all = init_rms(params.s_init)
    rms_loc = init_rms(params.s_init)

    #saving old data to see differences
    a_t_old = all_traj.traj[1].a_t
    weights_old = zeros(length(net_to_vec(θa)),2)

    for i_iter = 1:params.n_iter
        ĝ_a *= 0
        ĝ_v *= 0

        #execute the parallel_episodes
        res[i_iter], prec[i_iter] = parallel_episodes(params,θa,θv,all_traj,rms_all)

        #Reorganice Data in minibatches
        mini_batch = build_minibatches(all_traj,params)

        #First refit the value function for new data batch:
        θv = refit_VF(params, θv, mini_batch, opt_para_v)

        #Now apply the gae for the new data based on updated value function
        gae_all!(all_traj,θv,params)

        #Reorganice Data in minibatches as advantages have changed...
        mini_batch = build_minibatches(all_traj,params)

        #Do the Policy / and Value updates
        ĝ_a = parallel_gradient(params, θa, mini_batch)
        θa  = do_TNPG_update(ĝ_a, mini_batch, θa, params)

        if min_cost > res[i_iter]
            min_cost = res[i_iter]
            best_net = copy(θa)
            print(" new_best: ", round(min_cost, digits = 3))
        end

        #Visualization_stuff
        #-----------------------------------------------------------------------
        ProgressMeter.update!(progress,i_iter) #update progress bar
        if floor(i_iter/50) == ceil(i_iter/50)
            a_t_old ,weights_old = plot_stuff(params, all_traj, a_t_old,
            weights_old, i_iter,θa, θv, res, prec/params.n_traj, ĝ_a, ĝ_v)
        end
        #-----------------------------------------------------------------------

        @everywhere GC.gc(true)
    end
    return res
end

"""
TRPO
"""
function trpo(params::Params, θa::Any, θv::Any, res::Array{Float64,1}, prec::Array{Float64,1}, all_traj::Traj_Data, best_net::Any)
    #Init Part:
    best_net = copy(θa)
    min_cost = Inf
    #------------------------
    #Initialize the Policy
    ĝ_a = copy(θa)
    opt_para_a = optimizers(θa,   Adam; lr = params.lr_pl)

    #Initialize the Value Function
    ĝ_v = copy(θv)
    opt_para_v = optimizers(θv,   Adam; lr = params.lr_vf)

    #Init Progress bar
    progress = Progress(params.n_iter,1)

    # Init rms
    rms_all = init_rms(params.s_init)
    rms_loc = init_rms(params.s_init)

    #saving old data to see differences
    a_t_old = all_traj.traj[1].a_t
    weights_old = zeros(length(net_to_vec(θa)),2)

    for i_iter = 1:params.n_iter
        ĝ_a *= 0
        ĝ_v *= 0

        #execute the parallel_episodes
        res[i_iter], prec[i_iter] = parallel_episodes(params,θa,θv,all_traj,rms_all)

        #Reorganice Data in minibatches
        mini_batch = build_minibatches(all_traj,params)

        #First refit the value function for new data batch:
        θv = refit_VF(params, θv, mini_batch, opt_para_v)

        #Now apply the gae for the new data based on updated value function
        gae_all!(all_traj,θv,params)

        #Reorganice Data in minibatches as advantages have changed...
        mini_batch = build_minibatches(all_traj,params)

        #Apply the gradient
        ĝ_a = parallel_gradient(params, θa, mini_batch)
        θa  = do_TRPO_update(ĝ_a, mini_batch, θa, params)

        if min_cost > res[i_iter]
            min_cost = res[i_iter]
            best_net = copy(θa)
            print(" new_best: ", round(min_cost, digits = 3))
        end

        #Visualization_stuff
        #-----------------------------------------------------------------------
        ProgressMeter.update!(progress,i_iter) #update progress bar
        if floor(i_iter/50) == ceil(i_iter/50)
            a_t_old ,weights_old = plot_stuff(params, all_traj, a_t_old,
            weights_old, i_iter,θa, θv, res, prec/params.n_traj, ĝ_a, ĝ_v)
        end
        #-----------------------------------------------------------------------

        @everywhere GC.gc(true)
    end
    return res
end

function ppo(params::Params, θa::Any, θv::Any, res::Array{Float64,1}, prec::Array{Float64,1}, all_traj::Traj_Data, best_net::Any, i_start)
    #Init Part:
    best_net = copy(θa)
    min_cost = Inf
    #------------------------
    #Initialize the Policy
    ĝ_a = copy(θa)
    opt_para_a = optimizers(θa,   Adam; lr = params.lr_pl)

    #Initialize the Value Function
    ĝ_v = copy(θv)
    opt_para_v = optimizers(θv,   Adam; lr = params.lr_vf)

    #Init Progress bar
    progress = Progress(params.n_iter,1)

    # Init rms
    rms_all = init_rms(params.s_init)
    rms_loc = init_rms(params.s_init)

    a_t_old = all_traj.traj[1].a_t
    weights_old = zeros(length(net_to_vec(θa)),2)

    kl_means = zeros(params.n_iter)

    for i_iter = i_start:1:params.n_iter
        ĝ_a *= 0
        ĝ_v *= 0

        #Save data in intervalls
        if floor(i_iter/50) == ceil(i_iter/50)
            #save_data(res, prec, θa, θv, all_traj, i_iter)
        end

        #-----------------------------------------------------------------------
        #Solve it in Parallel!
        res[i_iter], prec[i_iter] = parallel_episodes(params,θa,θv,all_traj,rms_all)

        #Reorganice Data in minibatches
        mini_batch = build_minibatches(all_traj,params);

        #First refit the value function for new data batch:
        θv = refit_VF(params, θv, mini_batch, opt_para_v);

        #Now apply the gae for the new data based on updated value function
        gae_all!(all_traj,θv,params);

        #Reorganice Data in minibatches as advantages have changed...
        mini_batch = build_minibatches(all_traj,params)

        #Do the Policy / and Value updates
        (θa, dump) = parallel_PPO(params, θa, θv, mini_batch, opt_para_a, opt_para_v)

        #Get KL-div of last update step
        kl_means[i_iter] = kl_all(θa,all_traj::Traj_Data,params::Params)

        if min_cost > res[i_iter]
            min_cost = res[i_iter]
            best_net = copy(θa)
            print(" new_best: ", round(min_cost, digits = 3))
        end

        #Visualization_stuff
        #-----------------------------------------------------------------------
        ProgressMeter.update!(progress,i_iter) #update progress bar
        if floor(i_iter/50) == ceil(i_iter/50)
            a_t_old ,weights_old = plot_stuff(params, all_traj, a_t_old,
            weights_old, i_iter,θa, θv, res, prec/params.n_traj, ĝ_a, ĝ_v)
        end
        #-----------------------------------------------------------------------

        #clean the workspace from junk
        @everywhere GC.gc(true)
    end
    return res, kl_means
end

function ppo2(params::Params, θa::Any, θv::Any, res::Array{Float64,1}, prec::Array{Float64,1}, all_traj::Traj_Data, best_net::Any, i_start)
    #Init Part:
    best_net = copy(θa)
    min_cost = Inf
    #------------------------
    #Initialize the Policy
    ĝ_a = copy(θa)
    opt_para_a = optimizers(θa,   Adam; lr = params.lr_pl)

    #Initialize the Value Function
    ĝ_v = copy(θv)
    opt_para_v = optimizers(θv,   Adam; lr = params.lr_vf)

    #Init Progress bar
    progress = Progress(params.n_iter,1)

    # Init rms
    rms_all = init_rms(params.s_init)
    rms_loc = init_rms(params.s_init)

    a_t_old = all_traj.traj[1].a_t
    weights_old = zeros(length(net_to_vec(θa)),2)

    kl_means = zeros(params.n_iter)

    for i_iter = i_start:params.n_iter
        ĝ_a *= 0
        ĝ_v *= 0

        #Save data in intervalls
        if floor(i_iter/50) == ceil(i_iter/50)
            #save_data(res, prec, θa, θv, all_traj, i_iter)
        end

        #-----------------------------------------------------------------------
        #Solve it in Parallel!
        res[i_iter], prec[i_iter] = parallel_episodes(params,θa,θv,all_traj,rms_all)

        #Reorganice Data in minibatches
        mini_batch = build_minibatches(all_traj,params);

        #First refit the value function for new data batch:
        θv = refit_VF(params, θv, mini_batch, opt_para_v);

        #Now apply the gae for the new data based on updated value function
        gae_all!(all_traj,θv,params);

        #Reorganice Data in minibatches as advantages have changed...
        mini_batch = build_minibatches(all_traj,params)

        #Do the Policy / and Value updates
        (θa, dump) = parallel_PPO2(params, θa, θv, mini_batch, opt_para_a, opt_para_v)

        #Get KL-div of last update step
        kl_means[i_iter] = kl_all(θa,all_traj::Traj_Data,params::Params)

        if min_cost > res[i_iter]
            min_cost = res[i_iter]
            best_net = copy(θa)
            print(" new_best: ", round(min_cost, digits = 3))
        end

        #Visualization_stuff
        #-----------------------------------------------------------------------
        ProgressMeter.update!(progress,i_iter) #update progress bar
        if floor(i_iter/50) == ceil(i_iter/50)
            a_t_old ,weights_old = plot_stuff(params, all_traj, a_t_old,
            weights_old, i_iter,θa, θv, res, prec/params.n_traj, ĝ_a, ĝ_v)
        end
        #-----------------------------------------------------------------------

        #clean the workspace from junk
        @everywhere GC.gc(true)
    end
    return res, kl_means
end

"""
Time Delayed Deep Deterministic Policy Gradient
"""

include("src/td3/episodes.jl")
include("src/td3/gradients.jl")
include("src/td3/updates.jl")

#1. Randomly Initialize
#-------------------------------------------------------------

#1.3 Policy network
Pl           = init_θ_single_norm(params.structure_net_pl)
Pl_tar       = copy(Pl)
opt_para_pl  = optimizers(Pl, Adam; lr = params.lr_pl)

#1.4 Replay Buffer and Latest Trajectories
buffer           = Traj_Data()
buffer.traj      = Array{Trajectory,1}(undef,0)
buffer.timesteps = 0
buffer.n_traj    = 0

#1.5 Buffer for the most recent trajectories
all_traj  = init_all_traj(params)

#main function running DDPG
function td3(Pl::Any, res)

    #1.1 First Critic
    Qf_1         = init_θ_single_norm(params.structure_net_qf)
    Qf_tar_1     = copy(Qf_1)
    opt_para_qf1 = optimizers(Qf_1, Adam; lr = 0.5*params.lr_qf)

    #1.1 Second Critic
    Qf_2         = init_θ_single_norm(params.structure_net_qf2)
    Qf_tar_2     = copy(Qf_2)
    opt_para_qf2 = optimizers(Qf_2, Adam; lr = 2.0*params.lr_qf)

    #1.3 Thrid Actor
    #Pl           = init_θ_single_norm(params.structure_net_pl)
    Pl_tar       = copy(Pl)
    opt_para_pl  = optimizers(Pl, Adam; lr = params.lr_pl)

    weights_old = zeros(length(net_to_vec(Pl)),2)
    a_t_old = all_traj.traj[1].a_t

    #1.4 Init running mean std
    rms_all = init_rms(params.s_init)

    #Iterate
    for i_iter = 1:params.n_iter
        #Init Progress bar
        progress = Progress(params.n_iter,1)

        if floor(i_iter/100) == ceil(i_iter/100)
            #save_data(res, prec, Pl, Qf_1, all_traj, i_iter)
        end

        #2.Collect Data using current stochastic Policy
        #-------------------------------------------------------------
        #res[i_iter] = parallel_episodes_td3(params, Pl, Qf_1, Qf_2, all_traj,rms_all)
        res[i_iter], dump = parallel_episodes(params,Pl,Pl_tar,all_traj,rms_all)

        #3.Pump data in Replay Buffer:
        #-------------------------------------------------------------
        #3.1 Generate minibatches from current data_batch
        mini_batch = build_minibatches(all_traj,params)
        #3.2 Pump the data into the buffer
        pump_buffer!(buffer,mini_batch,params)

        #4.Train all the stuff now:
        #-------------------------------------------------------------
        do_td3_update_parallel!(Qf_1, Qf_tar_1, Qf_2, Qf_tar_2, Pl, Pl_tar, buffer, params, opt_para_pl, opt_para_qf1, opt_para_qf2)

        #5.Visualize Progress
        ProgressMeter.update!(progress,i_iter) #update progress bar

        #6. Plot data of interest
        if floor(i_iter/50) == ceil(i_iter/50)
            print("  new visualization: ")
            a_t_old ,weights_old = plot_stuff(params, all_traj, a_t_old,
            weights_old, i_iter,Pl, Qf_1, res, prec/params.n_traj, Pl, Pl)
        end
        print("  Result:  ", round(res[i_iter], digits = 3))
    end
    return res
end

#"""
#-==============================================================================
#                                 Get running
#-==============================================================================
#"""
using FileIO

testruns      = 1
#params.n_iter = 20

res      = zeros(params.n_iter)
prec     = zeros(params.n_iter)
all_traj = init_all_traj(params)
best_net = init_θ_single_norm(params.structure_net_pl)

filename = "ballbot"
openfile = string("Results/10_20_init_nets/",filename,".jld2")
init_nets = load(openfile)
result = zeros(6,params.n_iter)
zw_save = "zw_save.jld2"

RESULT = []
result = zeros(6,params.n_iter)
function evaluate_all_methods()

    for i_all = 1: testruns
        result = zeros(6,params.n_iter)
        kl_res = zeros(6,params.n_iter)
        fig_label  = string("All Methods ", i_all)
        fig_label2 = string("All Methods kl ", i_all)
        figure(fig_label)
        clf()

        println("REINFORCE")
        init_nets1 = load(openfile)
        θa = init_θ_single_norm(params.structure_net_pl)
        θv = init_θ_single_norm(params.structure_net_vf)
        θa  = copy(init_nets1["θa"])
        θv  = copy(init_nets1["θv"])
        result[1,:] = reinforce(params, θa, θv, zeros(params.n_iter), prec, all_traj, best_net)
        figure(fig_label)
        plot(result[1,:], label="REINFORCE")
        newfile = string("Results/REINFORCE_",i_all,"_data.jld2")
        @save newfile θa, θv, all_traj, best_net, result

        println("A2C")
        init_nets2 = load(openfile)
        θa = init_θ_single_norm(params.structure_net_pl)
        θv = init_θ_single_norm(params.structure_net_vf)
        θa  = copy(init_nets2["θa"])
        θv  = copy(init_nets2["θv"])
        result[2,:], kl_res[2,:] = advantage_actor_critic(params, θa, θv, zeros(params.n_iter), prec, all_traj, best_net)
        figure(fig_label)
        plot(result[2,:], label="A2C")
        legend()
        newfile = string("Results/A2C_",i_all,"_data.jld2")
        @save newfile θa, θv, all_traj, best_net, result, kl_res

        println("TNPG")
        init_nets3 = load(openfile)
        θa = init_θ_single_norm(params.structure_net_pl)
        θv = init_θ_single_norm(params.structure_net_vf)
        θa  = copy(init_nets3["θa"])
        θv  = copy(init_nets3["θv"])
        result[3,:] = tnpg(params, θa, θv, zeros(params.n_iter), prec, all_traj, best_net)
        figure(fig_label)
        plot(result[3,:], label="TNPG")
        legend()
        newfile = string("Results/TNPG_",i_all,"_data.jld2")
        @save newfile θa, θv, all_traj, best_net, result

        println("TRPO")
        init_nets4 = load(openfile)
        θa = init_θ_single_norm(params.structure_net_pl)
        θv = init_θ_single_norm(params.structure_net_vf)
        θa  = copy(init_nets4["θa"])
        θv  = copy(init_nets4["θv"])
        result[4,:] = trpo(params, θa, θv, zeros(params.n_iter), prec, all_traj, best_net)
        figure(fig_label)
        plot(result[4,:], label="TRPO")
        legend()
        newfile = string("Results/TRPO_",i_all,"_data.jld2")
        @save newfile θa, θv, all_traj, best_net, result

        println("PPO-Hard")
        init_nets5 = load(openfile)
        θa = init_θ_single_norm(params.structure_net_pl)
        θv = init_θ_single_norm(params.structure_net_vf)
        θa  = copy(init_nets5[1])
        θv  = copy(init_nets5[2])
        result[5,:], kl_res[5,:] = ppo(params, θa, θv, res5 ,prec, all_traj, best_net, 1024)
        figure(fig_label)
        plot(result[5,:], label="PPO-Hard")
        legend()
        figure(fig_label2)
        plot(kl_res[5,:], label="PPO-Hard")
        legend()
        newfile = string("Results/PPO1_",i_all,"_data.jld2")
        @save newfile θa, θv, all_traj, best_net, result, kl_res

        println("PPO-Soft")
        init_nets6 = load(openfile)
        θa = init_θ_single_norm(params.structure_net_pl)
        θv = init_θ_single_norm(params.structure_net_vf)
        θa  = copy(init_nets5[1])
        θv  = copy(init_nets5[2])
        result[6,:], kl_res[6,:] = ppo2(params, θa, θv, res6, prec, all_traj, best_net, 1024)
        figure(fig_label)
        plot(result[6,:], label="PPO-Soft")
        legend()
        newfile = string("Results/PPO2_",i_all,"_data.jld2")
        @save newfile θa, θv, all_traj, best_net, result, kl_res

        savefile = string("Results/result",Dates.today(),"_iter_",i_all,"_zw.jld2")
        @save savefile result

        append!(RESULT,result)
    end
end

evaluate_all_methods()

savefile = string("Results/Actor_Critic_",Dates.today(),"done.jld2")
@save savefile RESULT

max_stars_repo_name: NikonPic/REINFORCEMENT-JULIA

max_issues_repo_name: NikonPic/REINFORCEMENT-JULIA

max_forks_repo_name: NikonPic/REINFORCEMENT-JULIA


text:
function prepare_LS(stp, x₀, d, τ₀, f₀, ∇f₀)
    # extract the nlp
    nlp = stp.pb
    # construct the line search model, which will be rebased at each iteration'current data
    ϕ = LSModel(nlp, x₀, d, τ₀, f₀, dot(∇f₀,d) )
    # instantiate stp, which will be adjusted at each iteration's current data
    ϕstp = LS_Stopping(ϕ, LSAtT(0.0),
                       optimality_check = (p,s) -> optim_check_LS(p,s), #Tanj: already set optim_check_LS(p,s,α,β)?
                       main_stp = stp,
                       max_iter = 40,
                       atol = 0.0, # to rely only on the Armijo-Wolfe conditions
                       rtol = 0.0, # otherwise, tolerance may prohibit convergence
                       unbounded_threshold = 1e100)

    return ϕ, ϕstp
end


text:
using Base.Test

test_handler(r::Test.Success) = println("Test success: $(r.expr)")
test_handler(r::Test.Failure) = println("Test failure: $(r.expr)")
test_handler(r::Test.Error) = println("Test error: $(r.expr)")

f() = throw(DomainError())
g() = throw(BoundsError())
h() = error("An error happened")

Test.with_handler(test_handler) do
    @test_throws DomainError f()
    @test_throws BoundsError g()
    @test 1 != 1
    @test 1 == 1
    @test h()
    @test_throws BoundsError f()
    @test_throws DomainError g()
end

max_stars_repo_name: sjkelly/JuliaPlayground

max_issues_repo_name: sjkelly/JuliaPlayground

max_forks_repo_name: sjkelly/JuliaPlayground


text:
# Use baremodule to shave off a few KB from the serialized `.ji` file
baremodule NetCDFF_jll
using Base
using Base: UUID
import JLLWrappers

JLLWrappers.@generate_main_file_header("NetCDFF")
JLLWrappers.@generate_main_file("NetCDFF", UUID("78e728a9-57fe-5d11-897c-5014b89e5f84"))
end  # module NetCDFF_jll