text:
@testset "Is zero" begin
    @test is_zero([0])
    @test is_zero([0, 0])
    @test is_zero([0, 0, 0])

    @test !is_zero([1])
    @test !is_zero([1, 0])
    @test !is_zero([1, 0, 0])

    @test !is_zero([1e-2, 0, 0])
    @test is_zero([1e-2, 0, 0]; atol=1e-4)
    @test !is_zero([1e-1, 0, 0]; atol=1e-4)
end

@testset "Is Parallel" begin
    @test are_parallel([1, 0], [1, 0])
    @test are_parallel([1, 0], [-1, 0])
    @test are_parallel([1, 0], [2, 0])

    @test are_parallel([1, 5], [1, 5])
    @test are_parallel([1, 5], [3, 15])
    @test are_parallel([1, 5], [-5, -25])

    @test are_parallel([1, 2, 3], [3, 6, 9])

    @test !are_parallel([1, 0], [1, 1])
    @test !are_parallel([1, 0], [1, 1e-2])
    @test are_parallel([1, 0], [1, 1e-2]; atol=1e-2)

    # The zero vector is parallel to all vectors.
    @test are_parallel([0, 0], [1, 0])
    @test are_parallel([0, 0], [5, 3])
    @test are_parallel([3, 4], [0, 0])
    @test are_parallel([3, 4, -21], [0, 0, 0])

    @test !are_parallel([0, 1e-2], [1, 0])
    @test are_parallel([0, 1e-2], [1, 0]; atol=1e-2)
    @test !are_parallel([1, 0], [0, 1e-2])
    @test are_parallel([1, 0], [0, 1e-2]; atol=1e-2)
end

@testset "Is Perpendicular" begin
    @test are_perpendicular([1, 0], [0, 1])
    @test are_perpendicular([0, 1], [-1, 0])
    @test !are_perpendicular([1, 0], [1, 0])
    @test are_perpendicular([50, 0], [0, 2])

    @test are_perpendicular([5, 2], [2, -5])
    @test are_perpendicular([1, 1, 0], [0, 0, 1])

    @test !are_perpendicular([1, 0], [1e-2, 1])
    @test are_perpendicular([1, 0], [1e-2, 1]; atol=1e-2)
end

@testset "Are coplanar" begin

    # Any three or fewer points are coplanar.

    points = hcat([0; 0; 1])
    @test are_coplanar(points)

    points = [
        0 1;
        0 1;
        1 2;
    ]
    @test are_coplanar(points)

    points = [
        0 1 5;
        0 1 2;
        1 2 7;
    ]
    @test are_coplanar(points)

    # points = [
    #     0 0 1;
    #     1 1 0;
    #     5 2 0;
    #     6 -2 0;
    # ]
    points = [
        0 1 5 6;
        0 1 2 -2;
        1 0 0 0;
    ]
    @test !are_coplanar(points)

    points = [
        0 1 2 3;
        0 1 2 3;
        0 1 2 3;
    ]
    @test are_coplanar(points)

    points = [
        5 1 5 6;
        7 1 2 -2;
        4 4 4 4;
    ]
    @test are_coplanar(points)

    # Duplicates do not matter.
    points = [
        0 0 0;
        1 1 1;
        2 2 2;
        3 3 3;
        3 3 3;
    ]
    points = [
        0 1 2 3;
        0 1 2 3;
        0 1 2 3;
    ]
    @test are_coplanar(points)
end

@testset "Point on surface of line" begin
    for (point, line, bool_expected) in [
        ([0, 0], Line([0, 0], [1, 0]), true)
        ([2, 0], Line([0, 0], [5, 0]), true)
        ([2, 0], Line([0, 0], [-5, 0]), true)
        ([0, 0], Line([0, 0], [1, 1]), true)
        ([1, 0], Line([0, 0], [1, 1]), false)
        ([1, 0], Line([0, 0], [1, 1]), false)
    ]
        @test on_surface(point, line) == bool_expected
    end

    for (point, line, atol, bool_expected) in [
        ([0, 1e-2], Line([0, 0], [1, 0]), 1e-2, true)
        ([0, 1e-2], Line([0, 0], [1, 0]), 1e-3, false)
    ]
        @test on_surface(point, line; atol=atol) == bool_expected
    end
end

@testset "Lines are approximately equal" begin
    for (line_a, line_b, bool_expected) in [
        (Line([0, 0], [1, 0]), Line([0, 0], [1, 0]), true),
        (Line([0, 0], [1, 0]), Line([0, 0], [5, 0]), true),
        (Line([0, 0], [1, 0]), Line([0, 0], [-1, 0]), true),
        (Line([0, 0], [1, 0]), Line([0, 0], [-5, 0]), true),
        (Line([3, 0], [1, 0]), Line([0, 0], [-5, 0]), true),
        (Line([-3, 0], [1, 0]), Line([0, 0], [-5, 0]), true),
        (Line([0, 1], [1, 0]), Line([0, 0], [1, 0]), false),
        (Line([0, -1], [1, 0]), Line([0, 0], [1, 0]), false),
        (Line([0, 0], [1, 0]), Line([0, 0], [1, 1]), false),
    ]
        @test isapprox(line_a, line_b) == bool_expected
    end

    for (line_a, line_b, atol, bool_expected) in [
        (Line([0, 1e-2], [1, 0]), Line([0, 0], [1, 0]), 1e-2, true),
        (Line([0, 1e-2], [1, 0]), Line([0, 0], [1, 0]), 1e-3, false),
        (Line([0, 0], [1, 1]), Line([0, 0], [1, 1.01]), 1e-2, true),
        (Line([0, 0], [1, 1]), Line([0, 0], [1, 1.01]), 1e-5, false),
    ]
        @test isapprox(line_a, line_b; atol=atol) == bool_expected
    end
end


text:
eval_for(p::Poly, x₀) = p.a₀ + (p.n > 0 ? sum([coef*x₀^k for (k,coef) in p.terms]) : 0)
eval_for(p::Poly, x, x₀) = isequal(p.x, x) ? eval_for(p, x₀) : eval_for(sym(p), x, x₀)
eval_for(p, x, x₀) = substitute(p, Dict(x => x₀))

(p::Poly)(x₀) = eval_for(p, x₀)

"""
    solve_newton is a symbolic Newton-Ralphson solver
    f is a symbolic equation to be solved (f ~ 0)
    x is the variable to solve
    x₀ is the initial guess
"""
function solve_newton(f, x, x₀; abstol=1e-8, maxiter=50)
    xₙ = Complex(x₀)
    ∂f = value(derivative(f, x))

    for i = 1:maxiter
        xₙ₊₁ = xₙ - eval_for(f, x, xₙ) / eval_for(∂f, x, xₙ)

        if abs(xₙ₊₁ - xₙ) < abstol
            return xₙ₊₁
        else
            xₙ = xₙ₊₁
        end
    end
    return nothing
end

"""
    find_roots returns all the real and complex roots of a polynomial (poly)
    for variable x

    the output is rs, sz, where rs is a list of real roots and zs is a list of
    complex roots

    if poly is not a polynomial, the number of desired roots (n) needs to be
    provided. For example, find_roots(sin(x), x, 10) returns 10 different, but
    not necessarily sequential, multiples of π
"""
function find_roots(p, x, n=-1; abstol=1e-8)
    n = (n == -1 ? degree(p, x) : n)

    rs = Float64[]
    zs = Complex[]

    while n > 0
        z = solve_newton(p, x, Complex(rand(),rand()))
        if z != nothing
            if abs(imag(z)) < abstol
                r = real(z)
                push!(rs, r)
                if p isa Poly
                    p = p / gen_poly(x, -r, 1)
                else
                    p = p / (x - r)
                end
                n -= 1
            else
                if abs(real(z)) < abstol
                    z = Complex(0, imag(z))
                end
                append!(zs, [z, conj(z)])
                #poly = poly / ((x-z)*(x-conj(z)))
                if p isa Poly
                    p = p / gen_poly(x, abs2(z), -2*real(z), 1)
                else
                    p = p / (x^2 - 2x*real(z) + abs2(z))
                end
                n -= 2
            end
        end
    end
    sort(rs), zs
end

ext: jl


text:
#############################################################
#
# Streaming Functions
#
#############################################################

function post_streaming_filtered(; options=Dict{AbstractString, AbstractString}())
    error("Twitter API not fully implemented")
end

function get_streaming_sample(; options=Dict{AbstractString, AbstractString}())

    error("Twitter API not fully implemented")
    #r = get_oauth("https://stream.twitter.com/1.1/statuses/sample.json", options)

end

function get_firehose(; options=Dict{AbstractString, AbstractString}())
    error("This method requires special Twitter permission. If you have that special permission, write this function and submit a PR.")
end

function get_user_streaming(; options=Dict{AbstractString, AbstractString}())

    error("Twitter API not fully implemented")
    #r = get_oauth("https://userstream.twitter.com/1.1/user.json", options)

end

function get_site_streaming(; options=Dict{AbstractString, AbstractString}())

    error("Twitter API not fully implemented")
    #r = get_oauth("https://sitestream.twitter.com/1.1/site.json", options)

end


text:
"""
    set_properties!(dict, obj; defer_preview, opts...)

Internal function to set properties of an object `obj` given options (`opts`) of the form
`optname=value` an applying it through appropriate application function stored in the dictionary
`dict`. The keyword `defer_preview` is passed by functions that will apply a preview themselves
so that it's not repeated unnecessarily.
Note that the dictionary of property-setting-functions also contains "pre-conditioners" which
are functions that check that the values passed to properties are sensible and convert them
to sensible types if relevant. This reduces code duplication and allows to reduce specialization.
"""
function set_properties!(dict::Dict{Symbol,Pair{Function,Function}}, obj;
                         defer_preview=false, opts...)
    for optname ∈ opts.itr
        argcheck, setprop! = get(dict, optname) do
            throw(UnknownOptionError(optname, obj))
        end
        setprop!(obj, argcheck(opts[optname], optname))
    end
    defer_preview && return nothing
    return PreviewFigure(gcf())
end

set(obj; opts...) = set_properties!(obj; opts...)

####
#### Value checkers for set_properties functions the symbol corresponds to the name
#### of the option that is being modified
####

id(x, ::Symbol) = x
fl(x, ::Symbol) = fl(x)     # float conversion, see /utils.jl
not(x::Bool, ::Symbol) = !x # see for instance legend:nobox
lc(x::String, ::Symbol) = lowercase(x)
lc(x::Vector{String}, ::Symbol) = lowercase.(x)

"""
    posfl(x, s)

Internal function to check that all `x` are positive and then cast to Float64.
"""
function posfl(x, optname::Symbol)
    all(0 .< x) || throw(OptionValueError(String(optname), x))
    return fl(x)
end

"""
    posint(x, s)

Internal function to check that `x` is a positive integer.
"""
function posint(x::Int, optname::Symbol)
    0 < x || throw(OptionValueError(String(optname), x))
    return x
end

"""
    col(x, s)

Internal function to process a color, if a string is passed, it is parsed by the `Colors` package.
"""
col(c::Color, ::Symbol)   = c
col(s::String, ::Symbol)  = parse(Color, s)
col(v::Vector, s::Symbol) = col.(v, s)

"""
    opcol(x, s)

Internal function to process an optional color, i.e. `"none"` can be passed resulting in `nothing`
being passed. This is relevant for instance when setting the background color of a figure to "none"
(transparent figure). If `"none"` is passed, the current figure is checked for transparency
properties, if set to false, a warning will be raised and the color will be defaulted to `"white"`
if the transparency setting is unset, it will be set to `true`.
"""
function opcol(s::String, n::Symbol)::Option{Color}
    if lowercase(s) == "none"
        isdef(gcf().transparency) || (gcf().transparency=true)
        if !gcf().transparency
            @warn "Transparent background is only supported when the figure " *
                  "has its transparency property set to 'true'."
            return col("white", n) # fully opaque
        end
        return nothing
    end
    return col(s, n)
end
opcol(c::Color, n::Symbol) = col(c, n)

"""
    alpha(α, s)

Internal function to process an alpha parameter (transparency). The current figure is checked for
transparency properties, if set to false, a warning will be raised and the `α` will be ignored
(fully opaque). If the transparency setting is unset, it will be set to `true`.
Accepted values are strictly between `0` and `1`. For completely transparent, see using `"none"` in
the color description (for instance `Figure(bgcol="none")`).
"""
function alpha(α::Real, optname::Symbol)
    isdef(gcf().transparency) || (gcf().transparency=true)
    if !gcf().transparency
        @warn "Transparent colors are only supported when the figure " *
              "has its transparency property set to 'true'. Ignoring α."
        return 1.0 # fully opaque
    end
    0 < α < 1 || throw(OptionValueError("alpha"), α)
    return fl(α)
end

####
#### Pickers for GROUPED objects
#### --> gline, gbar, scatter, boxplot, bar
####

# this behaves like ∘ except it splats the output of g, this is useful in set_properties
# when wanting to use a set function but it has to be applied on a subfield of the object,
# then we apply a "selector" first to pick that subfield before applying the setter.
# the effect would be something like  setter(picker(obj), value) and is written (setter ⊙ picker)
⊙(f::Function, g::Function) = (λ(args...) = f(g(args...)...))

pick_lstyles(o, v) = (o.linestyles, v)
pick_mstyles(o, v) = (o.markerstyles, v)
pick_bstyles(o, v) = (o.barstyles, v)
pick_blstyle(b, v) = (b.boxstyles, v, :blstyle)
pick_mlstyle(b, v) = (b.boxstyles, v, :mlstyle)
pick_mmstyle(b, v) = (b.boxstyles, v, :mmstyle)
pick_omstyle(b, v) = (b.boxstyles, v, :omstyle)

###############################################################
####
#### Options for STYLE
####
###############################################################

const TEXTSTYLE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :font      => lc    => set_font!,      # set_style
    :fontsize  => posfl => set_hei!,       # .
    :col       => col   => set_textcolor!, # .
    :color     => col   => set_textcolor!, # .
    :textcol   => col   => set_textcolor!, # .
    :textcolor => col   => set_textcolor!, # .
    )

const LINESTYLE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :ls        => id    => set_lstyle!, # set_style
    :lstyle    => id    => set_lstyle!, # .
    :linestyle => id    => set_lstyle!, # .
    :lw        => posfl => set_lwidth!, # .
    :lwidth    => posfl => set_lwidth!, # .
    :linewidth => posfl => set_lwidth!, # .
    :smooth    => id    => set_smooth!, # .
    :col       => col   => set_color!,  # .
    :color     => col   => set_color!,  # .
    )

const GLINESTYLE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :ls         => id    => set_lstyles! ⊙ pick_lstyles, # set_style
    :lstyle     => id    => set_lstyles! ⊙ pick_lstyles, # .
    :linestyle  => id    => set_lstyles! ⊙ pick_lstyles, # .
    :lstyles    => id    => set_lstyles! ⊙ pick_lstyles, # .
    :linestyles => id    => set_lstyles! ⊙ pick_lstyles, # .
    :lw         => posfl => set_lwidths! ⊙ pick_lstyles, # .
    :lwidth     => posfl => set_lwidths! ⊙ pick_lstyles, # .
    :linewidth  => posfl => set_lwidths! ⊙ pick_lstyles, # .
    :lwidths    => posfl => set_lwidths! ⊙ pick_lstyles, # .
    :linewidths => posfl => set_lwidths! ⊙ pick_lstyles, # .
    :smooth     => id    => set_smooths! ⊙ pick_lstyles, # .
    :smooths    => id    => set_smooths! ⊙ pick_lstyles, # .
    :lc         => col   => set_colors! ⊙ pick_lstyles,  # .
    :col        => col   => set_colors! ⊙ pick_lstyles,  # .
    :color      => col   => set_colors! ⊙ pick_lstyles,  # .
    :cols       => col   => set_colors! ⊙ pick_lstyles,  # .
    :colors     => col   => set_colors! ⊙ pick_lstyles,  # .
    )

const MARKERSTYLE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :marker           => lc    => set_marker!, # set_style
    :ms               => posfl => set_marker!, # .
    :msize            => posfl => set_msize!,  # .
    :markersize       => posfl => set_msize!,  # .
    :mc               => col   => set_mcol!,   # .
    :mcol             => col   => set_mcol!,   # .
    :markercol        => col   => set_mcol!,   # .
    :markercolor      => col   => set_mcol!,   # .
    )

const GMARKERSTYLE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :marker           => lc    => set_markers! ⊙ pick_mstyles, # set_style
    :markers          => lc    => set_markers! ⊙ pick_mstyles, # .
    :ms               => posfl => set_markers! ⊙ pick_mstyles, # .
    :msize            => posfl => set_msizes! ⊙ pick_mstyles,  # .
    :msizes           => posfl => set_msizes! ⊙ pick_mstyles,  # .
    :markersize       => posfl => set_msizes! ⊙ pick_mstyles,  # .
    :markersizes      => posfl => set_msizes! ⊙ pick_mstyles,  # .
    :mc               => col   => set_mcols! ⊙ pick_mstyles,   # .
    :mcol             => col   => set_mcols! ⊙ pick_mstyles,   # .
    :markercol        => col   => set_mcols! ⊙ pick_mstyles,   # .
    :markercolor      => col   => set_mcols! ⊙ pick_mstyles,   # .
    :mcols            => col   => set_mcols! ⊙ pick_mstyles,   # .
    :markercols       => col   => set_mcols! ⊙ pick_mstyles,   # .
    :markercolors     => col   => set_mcols! ⊙ pick_mstyles,   # .
    )

const BARSTYLE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :ecol      => col => set_color!, # .
    :edgecol   => col => set_color!, # .
    :edgecolor => col => set_color!, # .
    )

const GBARSTYLE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    # bar fill color (main)
    :col        => col   => set_fills! ⊙ pick_bstyles,  # set_style
    :color      => col   => set_fills! ⊙ pick_bstyles,  # .
    :cols       => col   => set_fills! ⊙ pick_bstyles,  # .
    :colors     => col   => set_fills! ⊙ pick_bstyles,  # .
    :fill       => col   => set_fills! ⊙ pick_bstyles,  # .
    :fills      => col   => set_fills! ⊙ pick_bstyles,  # .
    # edge color (secondary)
    :ecol       => col   => set_colors! ⊙ pick_bstyles, # .
    :edgecol    => col   => set_colors! ⊙ pick_bstyles, # .
    :edgecolor  => col   => set_colors! ⊙ pick_bstyles, # .
    :ecols      => col   => set_colors! ⊙ pick_bstyles, # .
    :edgecols   => col   => set_colors! ⊙ pick_bstyles, # .
    :edgecolors => col   => set_colors! ⊙ pick_bstyles, # .
    # bar width
    :width      => posfl => set_bwidth!, # .
    :bwidth     => posfl => set_bwidth!, # .
    :barwidth   => posfl => set_bwidth!, # .
    )

const FILLSTYLE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    # fill color and alpha
    :col   => col   => set_fill!,  # set_style
    :color => col   => set_fill!,  # .
    :fill  => col   => set_fill!,  # .
    :alpha => alpha => set_alpha!, # .
    )

###############################################################
####
#### Options for AX_ELEMS
####
###############################################################

const TITLE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :dist => posfl => set_dist!
    )
merge!(TITLE_OPTS, TEXTSTYLE_OPTS)
set_properties!(t::Title; opts...) = set_properties!(TITLE_OPTS, t; opts...)

const LEGEND_OPTS = Dict{Symbol,Pair{Function,Function}}(
    # position of the legend
    :pos        => lc    => set_position!, # set_legend
    :position   => lc    => set_position!, # .
    # toggle legend on/off
    :off        => id    => set_off!,      # .
    :on         => not   => set_off!,      # .
    # toggle surrounding box on/off
    :nobox      => id    => set_nobox!,    # .
    :box        => not   => set_nobox!,    # .
    # specify margins (internal) and offset (external)
    :margins    => fl    => set_margins!,  # .
    :offset     => fl    => set_offset!,   # .
    # background color and alpha
    :bgcol      => opcol => set_color!,    # set_style
    :bgcolor    => opcol => set_color!,    # .
    :background => opcol => set_color!,    # .
    :bgalpha    => alpha => set_alpha!,    # .
    )
merge!(LEGEND_OPTS, TEXTSTYLE_OPTS)
set_properties!(l::Legend; opts...) = set_properties!(LEGEND_OPTS, l; opts...)

const TICKS_OPTS = Dict{Symbol,Pair{Function,Function}}(
    # ticks related
    :off        => id  => set_off!,        # set_ax_elems
    :on         => not => set_off!,        # .
    :len        => fl  => set_length!,     # .
    :length     => fl  => set_length!,     # .
    :sym        => id  => set_symticks!,   # .
    :symticks   => id  => set_symticks!,   # .
    # ticks color
    :tcol       => col => set_color!,      # .
    :tickscol   => col => set_color!,      # .
    :tickscolor => col => set_color!,      # .
    # grid mode
    :grid       => id  => set_grid!,       # .
    # labels related
    :hidelabels => id  => set_labels_off!, # set_ax_elems
    :showlabels => not => set_labels_off!, # .
    :angle      => fl  => set_angle!,      # .
    :format     => lc  => set_format!,     # .
    :shift      => fl  => set_shift!,      # .
    :dist       => id  => set_dist!,       # .
    :distance   => id  => set_dist!,
    )
merge!(TICKS_OPTS, LINESTYLE_OPTS) # ticks line
merge!(TICKS_OPTS, TEXTSTYLE_OPTS) # labels
set_properties!(t::Ticks; opts...) = set_properties!(TICKS_OPTS, t; opts...)

###############################################################
####
#### Options for DRAWINGS
####
###############################################################

const SCATTER2D_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :key    => id => set_labels!, # .
    :keys   => id => set_labels!, # .
    :label  => id => set_labels!, # .
    :labels => id => set_labels!, # .
    )
merge!(SCATTER2D_OPTS, GLINESTYLE_OPTS)
merge!(SCATTER2D_OPTS, GMARKERSTYLE_OPTS)
set_properties!(s::Scatter2D; opts...) = set_properties!(SCATTER2D_OPTS, s; opts...)

const FILL2D_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :from  => fl => set_xmin!,
    :min   => fl => set_xmin!,
    :xmin  => fl => set_xmin!,
    :to    => fl => set_xmax!,
    :max   => fl => set_xmax!,
    :xmax  => fl => set_xmax!,
    :key   => id => set_label!, # set_drawing
    :label => id => set_label!, # .
    )
merge!(FILL2D_OPTS, FILLSTYLE_OPTS)
set_properties!(f::Fill2D; opts...) = set_properties!(FILL2D_OPTS, f; opts...)

const HIST2D_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :bins       => posint => set_bins!,    # set_drawing
    :nbins      => posint => set_bins!,    # .
    :scaling    => lc     => set_scaling!, # .
    :norm       => lc     => set_scaling!, # .
    :horiz      => id     => set_horiz!,   # .
    :horizontal => id     => set_horiz!,   # .
    :key        => id     => set_label!,   # set_drawing
    :label      => id     => set_label!,   # .
    )
merge!(HIST2D_OPTS, BARSTYLE_OPTS)
merge!(HIST2D_OPTS, FILLSTYLE_OPTS)
set_properties!(h::Hist2D; opts...) = set_properties!(HIST2D_OPTS, h; opts...)

const BAR2D_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :stacked    => id => set_stacked!, # set_drawing
    :horiz      => id => set_horiz!,   # .
    :horizontal => id => set_horiz!,   # .
    :key        => id => set_labels!,  # .
    :label      => id => set_labels!,  # .
    :keys       => id => set_labels!,  # .
    :labels     => id => set_labels!,  # .
    )
merge!(BAR2D_OPTS, GBARSTYLE_OPTS)
set_properties!(gb::Bar2D; opts...) = set_properties!(BAR2D_OPTS, gb; opts...)

const BOXPLOT_OPTS = Dict{Symbol,Pair{Function,Function}}(
    # -- global toggle
    :horiz      => id => set_horiz!, # set_drawing
    :horizontal => id => set_horiz!,
    # -- now inside BOXSTYLES
    # box styling
    :box_width          => posfl => set_bwidths!, # set_drawing
    :box_widths         => posfl => set_bwidths!, # .
    :box_wwidth         => posfl => set_wwidths!, # .
    :box_wwidths        => posfl => set_wwidths!, # .
    :box_whisker_width  => posfl => set_wwidths!, # .
    :box_whisker_widths => posfl => set_wwidths!, # .
    # how long should the whiskers be
    :whisker            => posfl => set_wrlengths!, # set_drawing
    :whiskers           => posfl => set_wrlengths!, # .
    # what line style should be used to draw the boxes
    :box_ls         => id    => set_lstyles! ⊙ pick_blstyle,
    :box_lstyle     => id    => set_lstyles! ⊙ pick_blstyle,
    :box_lstyles    => id    => set_lstyles! ⊙ pick_blstyle,
    :box_linestyle  => id    => set_lstyles! ⊙ pick_blstyle,
    :box_linestyles => id    => set_lstyles! ⊙ pick_blstyle,
    :box_lw         => posfl => set_lwidths! ⊙ pick_blstyle,
    :box_lwidth     => posfl => set_lwidths! ⊙ pick_blstyle,
    :box_lwidths    => posfl => set_lwidths! ⊙ pick_blstyle,
    :box_linewidth  => posfl => set_lwidths! ⊙ pick_blstyle,
    :box_linewidths => posfl => set_lwidths! ⊙ pick_blstyle,
    :box_col        => col   => set_colors! ⊙ pick_blstyle,
    :box_cols       => col   => set_colors! ⊙ pick_blstyle,
    :box_color      => col   => set_colors! ⊙ pick_blstyle,
    :box_colors     => col   => set_colors! ⊙ pick_blstyle,
    # median line
    :med_ls         => id    => set_lstyles! ⊙ pick_mlstyle,
    :med_lstyle     => id    => set_lstyles! ⊙ pick_mlstyle,
    :med_lstyles    => id    => set_lstyles! ⊙ pick_mlstyle,
    :med_linestyle  => id    => set_lstyles! ⊙ pick_mlstyle,
    :med_linestyles => id    => set_lstyles! ⊙ pick_mlstyle,
    :med_lw         => posfl => set_lwidths! ⊙ pick_mlstyle,
    :med_lwidth     => posfl => set_lwidths! ⊙ pick_mlstyle,
    :med_lwidths    => posfl => set_lwidths! ⊙ pick_mlstyle,
    :med_linewidth  => posfl => set_lwidths! ⊙ pick_mlstyle,
    :med_linewidths => posfl => set_lwidths! ⊙ pick_mlstyle,
    :med_col        => col   => set_colors! ⊙ pick_mlstyle,
    :med_cols       => col   => set_colors! ⊙ pick_mlstyle,
    :med_color      => col   => set_colors! ⊙ pick_mlstyle,
    :med_colors     => col   => set_colors! ⊙ pick_mlstyle,
    # mean marker
    :mean_show         => id    => set_mshow!,
    :mean_marker       => id    => set_markers! ⊙ pick_mmstyle,
    :mean_markers      => id    => set_markers! ⊙ pick_mmstyle,
    :mean_msize        => posfl => set_msizes! ⊙ pick_mmstyle,
    :mean_msizes       => posfl => set_msizes! ⊙ pick_mmstyle,
    :mean_markersize   => posfl => set_msizes! ⊙ pick_mmstyle,
    :mean_markersizes  => posfl => set_msizes! ⊙ pick_mmstyle,
    :mean_mcol         => col   => set_mcols! ⊙ pick_mmstyle,
    :mean_mcols        => col   => set_mcols! ⊙ pick_mmstyle,
    :mean_markercol    => col   => set_mcols! ⊙ pick_mmstyle,
    :mean_markercols   => col   => set_mcols! ⊙ pick_mmstyle,
    :mean_markercolor  => col   => set_mcols! ⊙ pick_mmstyle,
    :mean_markercolors => col   => set_mcols! ⊙ pick_mmstyle,
    # outliers
    :out_show         => id    => set_mshow!,
    :out_marker       => id    => set_markers! ⊙ pick_omstyle,
    :out_markers      => id    => set_markers! ⊙ pick_omstyle,
    :out_msize        => posfl => set_msizes! ⊙ pick_omstyle,
    :out_msizes       => posfl => set_msizes! ⊙ pick_omstyle,
    :out_markersize   => posfl => set_msizes! ⊙ pick_omstyle,
    :out_markersizes  => posfl => set_msizes! ⊙ pick_omstyle,
    :out_mcol         => col   => set_mcols! ⊙ pick_omstyle,
    :out_mcols        => col   => set_mcols! ⊙ pick_omstyle,
    :out_markercol    => col   => set_mcols! ⊙ pick_omstyle,
    :out_markercols   => col   => set_mcols! ⊙ pick_omstyle,
    :out_markercolor  => col   => set_mcols! ⊙ pick_omstyle,
    :out_markercolors => col   => set_mcols! ⊙ pick_omstyle,
    )
set_properties!(bp::Boxplot; opts...) = set_properties!(BOXPLOT_OPTS, bp; opts...)

const HEATMAP_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :cmap     => col => set_cmap!,
    :colormap => col => set_cmap!,
    )
set_properties!(h::Heatmap; opts...) = set_properties!(HEATMAP_OPTS, h; opts...)

const SCATTER3D_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :key    => id => set_label!, # .
    :label  => id => set_label!, # .
    )
merge!(SCATTER3D_OPTS, LINESTYLE_OPTS)
merge!(SCATTER3D_OPTS, MARKERSTYLE_OPTS)
set_properties!(s::Scatter3D; opts...) = set_properties!(SCATTER3D_OPTS, s; opts...)

###############################################################
####
#### Options for OBJECTS
####
###############################################################

const TEXT2D_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :pos      => id    => set_position!, # set_drawing
    :position => id    => set_position!, # .
    )
merge!(TEXT2D_OPTS, TEXTSTYLE_OPTS)
set_properties!(t::Text2D; opts...) = set_properties!(TEXT2D_OPTS, t; opts...)

const STRAIGHTLINE2D_OPTS = Dict{Symbol,Pair{Function,Function}}(
    )
merge!(STRAIGHTLINE2D_OPTS, LINESTYLE_OPTS)
set_properties!(t::StraightLine2D; opts...) = set_properties!(STRAIGHTLINE2D_OPTS, t; opts...)

const BOX2D_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :pos      => id    => set_position!, # set legend
    :position => id    => set_position!, # .
    :nobox    => id    => set_nobox!,    # .
    :box      => not   => set_nobox!,    # .
    :fill     => opcol => set_fill!,     # .
    :alpha    => alpha => set_alpha!,    # .
    )
merge!(BOX2D_OPTS, LINESTYLE_OPTS)
set_properties!(b::Box2D; opts...) = set_properties!(BOX2D_OPTS, b; opts...)

const COLORBAR_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :pixels     => posint => set_pixels!,   # set object
    :res        => posint => set_pixels!,   # .
    :resolution => posint => set_pixels!,   # .
    :nobox      => id     => set_nobox!,    # set legend
    :box        => not    => set_nobox!,    # .
    :pos        => id     => set_position!, # set object
    :position   => id     => set_position!, # .
    :offset     => fl     => set_offset!,   # set legend
    :size       => posfl  => set_size!,     # set figure
    :ticks      => fl     => set_ticks!,    # set object
    :labels     => id     => set_labels!,
    )
merge!(COLORBAR_OPTS, TICKS_OPTS)
merge!(COLORBAR_OPTS, TEXTSTYLE_OPTS)
set_properties!(b::Colorbar; opts...) = set_properties!(COLORBAR_OPTS, b; opts...)

###############################################################
####
#### Options for AX*
####
###############################################################

const  AXES_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :size  => fl => set_size!,  # set figure
    :title => id => set_title!, # set ax
    :off   => id => set_off!,   # set axelems
    )
set_properties!(a::Axes; opts...) = set_properties!(AXES_OPTS, a; opts...)

const AXIS_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :title  => id    => set_title!,  # set_ax
    :base   => posfl => set_base!,   # .
    :min    => fl    => set_min!,    # .
    :max    => fl    => set_max!,    # .
    :log    => id    => set_log!,    # .
    :lwidth => posfl => set_lwidth!, # set_style
    :off    => id    => set_off!,    # set_ax_elems
    )
merge!(AXIS_OPTS, TEXTSTYLE_OPTS)
set_properties!(a::Axis; opts...) = set_properties!(AXIS_OPTS, a; opts...)

###############################################################
####
#### Options for FIGURE
####
###############################################################

const FIGURE_OPTS = Dict{Symbol,Pair{Function,Function}}(
    :size         => posfl => set_size!,         # set_figure
    :tex          => id    => set_texlabels!,    # .
    :hastex       => id    => set_texlabels!,    # .
    :latex        => id    => set_texlabels!,    # .
    :haslatex     => id    => set_texlabels!,    # .
    :texscale     => lc    => set_texscale!,     # .
    :alpha        => id    => set_transparency!, # .
    :transparent  => id    => set_transparency!, # .
    :transparency => id    => set_transparency!, # .
    :preamble     => id    => set_texpreamble!,  # .
    :texpreamble  => id    => set_texpreamble!,  # .
    :col          => opcol => set_color!,        # set_style
    :color        => opcol => set_color!,        # .
    :bgcol        => opcol => set_color!,        # .
    :bgcolor      => opcol => set_color!,        # .
    :background   => opcol => set_color!,        # .
    :bgalpha      => alpha => set_alpha!,        # .
    )
merge!(FIGURE_OPTS, TEXTSTYLE_OPTS)
set_properties!(f::Figure; opts...) = set_properties!(FIGURE_OPTS, f; opts...)


text:
@testset "Hawkes Fit" begin

    testEvents = sort(rand(100)*100)

    testSamples = HawkesProcesses.fit(testEvents, 100, 3)
    inferTests && @inferred HawkesProcesses.fit(testEvents, 100, 3)
    @test isa(testSamples, Tuple)
    @test length(testSamples) == 3
    @test size(testSamples[1]) == (3,)
    @test size(testSamples[2]) == (3,)
    @test size(testSamples[3]) == (3,)

end


text:
using Documenter, TransformVariables

DocMeta.setdocmeta!(TransformVariables, :DocTestSetup, :(using TransformVariables); recursive=true)
makedocs(
    modules = [TransformVariables],
    format = Documenter.HTML(; prettyurls = get(ENV, "CI", nothing) == "true"),
    sitename = "TransformVariables.jl",
    clean = true,
    checkdocs = :exports,
    pages = Any["index.md", "internals.md"]
)

deploydocs(repo = "github.com/tpapp/TransformVariables.jl.git")


text:
# we write a function that returns the row index of the (absolutely) largest element in column k of the matrix A
function findLargestInCol(A, k)

	# we get column k
	col = A[:,k]

	# we take elementwise absolute value
	absCol = broadcast(abs, col)

	# we get the index of the largest element
	index = argmax(absCol)

	# we return index
	return index
end

max_stars_repo_path: 1a/LU with pivot/findLargestInCol.jl

max_issues_repo_path: 1a/LU with pivot/findLargestInCol.jl

max_forks_repo_path: 1a/LU with pivot/findLargestInCol.jl


text:
### Generate linear data through Born scattering modeling for each num_vintage
## author: Ziyi Yin (ziyi.yin@gatech.edu)

using DrWatson
@quickactivate "CompassTimeLapseCCS"

using PyPlot, JUDI
using JLD2, LinearAlgebra
using JOLI, FFTW
using Printf, Statistics

JLD2.@load "../model/CompassTimeLapse.jld2"
JLD2.@load "../model/MarineSourceLoc.jld2" xrec_stack zrec_stack
JLD2.@load "../model/OBN.jld2" OBN_loc

model0_stack = [Model(n,d,o,m0_stack[i]; nb=100) for i = 1:num_vintage]

### reciprocity

nsrc = length(OBN_loc)
nrec = length(xrec_stack[1])

dtS = dtR = 1f0
timeS = timeR = 2000f0

xsrc = convertToCell(OBN_loc)
ysrc = convertToCell(range(0f0,stop=0f0,length=nsrc))
zsrc = convertToCell(range(12*d[2],stop=12*d[2],length=nsrc))

srcGeometry = Geometry(xsrc, ysrc, zsrc; dt=dtS, t=timeS)
wavelet = ricker_wavelet(timeS,dtS,0.025f0)

q = judiVector(srcGeometry, wavelet)

yrec = 0f0

recGeometry_stack = [Geometry(xrec_stack[i],yrec,zrec_stack[i]; dt=dtR, t=timeR, nsrc=nsrc) for i = 1:num_vintage]

ntComp = get_computational_nt(srcGeometry, recGeometry_stack[1], model0_stack[1])
info = Info(prod(n), nsrc, ntComp)

opt = Options(isic=true)

F0_stack = [judiProjection(info, recGeometry_stack[i])*judiModeling(info, model0_stack[i]; options=opt)*adjoint(judiProjection(info, srcGeometry)) for i = 1:num_vintage]
J_stack = [judiJacobian(F0_stack[i],q) for i = 1:num_vintage]
d_lin = [J_stack[i]*dm_stack[i] for i = 1:num_vintage]    # generate linear data through Born scattering

JLD2.@save "../data/d_lin.jld2" d_lin q


text:
using FixedPointNumbers

@testset "Parse" begin
    r8(x) = reinterpret(N0f8, x)

    # Color parsing
    redN0f8 = parse(Colorant, "red")
    @test colorant"red" == redN0f8
    @test isa(redN0f8, RGB{N0f8})
    @test redN0f8 == RGB(1,0,0)
    @test parse(RGB{Float64}, "red") === RGB{Float64}(1,0,0)
    @test isa(parse(HSV, "blue"), HSV)
    @test parse(Colorant, "rgb(55,217,127)") === RGB{N0f8}(r8(0x37),r8(0xd9),r8(0x7f))
    @test colorant"rgb(55,217,127)" === RGB{N0f8}(r8(0x37),r8(0xd9),r8(0x7f))
    @test parse(Colorant, "rgba(55,217,127,0.5)") === RGBA{N0f8}(r8(0x37),r8(0xd9),r8(0x7f),0.5)
    @test parse(Colorant, "rgb(55,217,127)") === RGB{N0f8}(r8(0x37),r8(0xd9),r8(0x7f))
    @test parse(Colorant, "rgba(55,217,127,0.5)") === RGBA{N0f8}(r8(0x37),r8(0xd9),r8(0x7f),0.5)
    @test parse(Colorant, "hsl(120, 100%, 50%)") === HSL{Float32}(120,1.0,.5)
    @test colorant"hsl(120, 100%, 50%)" === HSL{Float32}(120,1.0,.5)
    @test parse(RGB{N0f8}, "hsl(120, 100%, 50%)") === convert(RGB{N0f8}, HSL{Float32}(120,1.0,.5))
    @test_throws ErrorException  parse(Colorant, "hsl(120, 100, 50)")
    @test parse(Colorant, "#D0FF58") === RGB(r8(0xD0),r8(0xFF),r8(0x58))
    @test parse(Colorant, "#FB0") === RGB(r8(0xFF),r8(0xBB),r8(0x00))

    @test parse(Colorant, :red) === colorant"red"
    @test parse(Colorant, colorant"red") === colorant"red"

end


text:
#####
##### three-dimensional fluxes
#####

"function euler_fluxes_3D(rhoL,uL,vL,wL,betaL,rhologL,betalogL,
                         rhoR,uR,vR,wR,betaR,rhologR,betalogR)"
function euler_fluxes_3D(rhoL,uL,vL,wL,betaL,rhologL,betalogL,
                         rhoR,uR,vR,wR,betaR,rhologR,betalogR)

    rholog = logmean.(rhoL,rhoR,rhologL,rhologR)
    betalog = logmean.(betaL,betaR,betalogL,betalogR)

    # arithmetic avgs
    rhoavg = (@. .5*(rhoL+rhoR))
    uavg   = (@. .5*(uL+uR))
    vavg   = (@. .5*(vL+vR))
    wavg   = (@. .5*(wL+wR))

    unorm = (@. uL*uR + vL*vR + wL*wR)
    pa    = (@. rhoavg/(betaL+betaR))
    E_plus_p  = (@. rholog/(2*(γ-1)*betalog) + pa + .5*rholog*unorm)

    FxS1 = (@. rholog*uavg)
    FxS2 = (@. FxS1*uavg + pa)
    FxS3 = (@. FxS1*vavg) # rho * u * v
    FxS4 = (@. FxS1*wavg) # rho * u * w
    FxS5 = (@. E_plus_p*uavg)

    FyS1 = (@. rholog*vavg)
    FyS2 = (@. FxS3) # rho * u * v
    FyS3 = (@. FyS1*vavg + pa)
    FyS4 = (@. FyS1*wavg) # rho * v * w
    FyS5 = (@. E_plus_p*vavg)

    FzS1 = (@. rholog*wavg)
    FzS2 = (@. FxS4) # rho * u * w
    FzS3 = (@. FyS4) # rho * v * w
    FzS4 = (@. FzS1*wavg + pa) # rho * w^2 + p
    FzS5 = (@. E_plus_p*wavg)

    Fx = (FxS1,FxS2,FxS3,FxS4,FxS5)
    Fy = (FyS1,FyS2,FyS3,FyS4,FyS5)
    Fz = (FzS1,FzS2,FzS3,FzS4,FzS5)
    return Fx,Fy,Fz
end

#####
##### shim functions
#####

# dispatch to n-dimensional constitutive routines, with optional entropy scaling
function primitive_to_conservative(rho,u,v,w,p)
   rho,rhoU,E = primitive_to_conservative_nd(rho,(u,v,w),p)
   return rho,rhoU...,E
end
function v_ufun(rho,rhou,rhov,rhow,E)
    v1,vU,vE = v_ufun_nd(rho,(rhou,rhov,rhow),E)
    return scale_entropy_output(v1,vU...,vE)
end
function u_vfun(v1,vU1,vU2,vU3,vE)
    v1,vU1,vU2,vU3,vE = scale_entropy_input(v1,vU1,vU2,vU3,vE)
    rho,rhoU,E = u_vfun_nd(v1,(vU1,vU2,vU3),vE)
    return rho,rhoU...,E
end
function conservative_to_primitive_beta(rho,rhou,rhov,rhow,E)
    rho,U,beta = conservative_to_primitive_beta_nd(rho,(rhou,rhov,rhow),E)
    return rho,U...,beta
end

Sfun(rho,rhou,rhov,rhow,E) = Sfun_nd(rho,(rhou,rhov,rhow),E)
pfun(rho,rhou,rhov,rhow,E) = pfun_nd(rho,(rhou,rhov,rhow),E)
betafun(rho,rhou,rhov,rhow,E) = betafun_nd(rho,(rhou,rhov,rhow),E)

"function euler_fluxes_3D(rhoL,uL,vL,wL,betaL,rhoR,uR,vR,wR,betaR)"
function euler_fluxes(rhoL,uL,vL,wL,betaL,rhoR,uR,vR,wR,betaR)
    rhologL,betalogL,rhologR,betalogR = map(x->log.(x),(rhoL,betaL,rhoR,betaR))
    return euler_fluxes_3D(rhoL,uL,vL,wL,betaL,rhologL,betalogL,
                           rhoR,uR,vR,wR,betaR,rhologR,betalogR)
end
euler_fluxes(rhoL,uL,vL,wL,betaL,rhologL,betalogL,
             rhoR,uR,vR,wR,betaR,rhologR,betalogR) =
     euler_fluxes_3D(rhoL,uL,vL,wL,betaL,rhologL,betalogL,
                     rhoR,uR,vR,wR,betaR,rhologR,betalogR)

max_stars_repo_name: yiminllin/ESDG-PosLimit

max_issues_repo_name: yiminllin/ESDG-PosLimit

max_forks_repo_name: yiminllin/ESDG-PosLimit


text:
using LinearAlgebra

abstract type AbstractMultiRobotProblem <: PartitionProblem{Tuple{Int64,
                                                                  Trajectory}}
end

function get_state(p::AbstractMultiRobotProblem, index)
  p.partition_matroid[index]
end

# Return agent center (which is just the state)
# The original coverage code involved more complex agent representations.
get_center(x::Tuple{Int64,Int64}) = [x...]

function sample_block(p::AbstractMultiRobotProblem, block::Integer)
  horizon = p.configs.horizon
  grid = p.configs.grid

  trajectory = Array{State}(undef, horizon)

  # (target dynamics are the same as random rollouts for tracking robots)
  current_state = p.partition_matroid[block]
  for ii in 1:horizon
    current_state = target_dynamics(grid, current_state)
    trajectory[ii] = current_state
  end

  (block, trajectory)
end

max_stars_repo_name: mcorah/MultiAgentSensing

max_issues_repo_name: mcorah/MultiAgentSensing

max_forks_repo_name: mcorah/MultiAgentSensing


text:
export
    BinaryOptimizationModel,
    vars,
    values,
    order,
    isquadratic

import Base: *, /, +, -, ^, getindex, setindex!

abstract type BinaryOptimizationModel{T<:Real} end

vars(bqm::BinaryOptimizationModel) = bqm.vars
values(bqm::BinaryOptimizationModel) = bqm.values

function _dict_add!(dict::Dict{Vector{Vector{Int}},T},
                    index::Vector{Vector{Int}},
                    val::Real) where T<:Real
    if index ∈ keys(dict)
        dict[index] += val
    else
        dict[index] = val
    end
    if dict[index] == zero(T)
        pop!(dict, index)
    end
    dict
end

function _dict_add_no_check!(dict::Dict{Vector{Vector{Int}},T},
                             index::Vector{Vector{Int}},
                             val::Real) where T<:Real
    if index ∈ keys(dict)
        dict[index] += val
    else
        dict[index] = val
    end
    dict
end

function getindex(bqm::S, key...) where S<:BinaryOptimizationModel{T} where T<:Real
    if isempty(key)
        key = Vector{Int}[]
    end
    the_key = clean_key(S, collect(key))
    if the_key ∈ keys(values(bqm))
        values(bqm)[the_key]
    else
        zero(T)
    end
end

function setindex!(bqm::S,
                   val::Real,
                   key...) where S<:BinaryOptimizationModel{T} where T<:Real
    if isempty(key)
        key = Vector{Int}[]
    end
    the_key = clean_key(S, collect(key))
    if val != zero(T)
        values(bqm)[the_key] = T(val)
    elseif the_key ∈ keys(values(bqm))
        delete!(values(bqm), the_key)
    end
    bqm
end

## Multiplication by number
function *(bqm::S, var::T) where {S<:BinaryOptimizationModel{T} where T<:Real, T<:Real}
    bqm_copy = S(vars(bqm), copy(values(bqm)))
    if var == zero(T)
        empty!(values(bqm_copy))
        return bqm_copy
    elseif var == one(T)
        return bqm_copy
    else
        for k = keys(values(bqm_copy))
            bqm_copy[k...] *= var
        end
    end
    bqm_copy
end

*(var::Real, bqm::BinaryOptimizationModel) = *(bqm, var)
/(bqm::BinaryOptimizationModel, var::T) where T<:Real = *(bqm, one(T)/var)
-(bqm::BinaryOptimizationModel{T}) where T = *(bqm, -one(T))
^(bqm::BinaryOptimizationModel, n::Int) = n == 1 ? bqm : bqm*bqm^(n-1)

## sum with number
function +(bqm::S, var::T) where {S<:BinaryOptimizationModel{T} where T<:Real,T<:Real}
    bqm_copy = S(vars(bqm), copy(values(bqm)))
    empty_array = Vector{Int}[]
    _dict_add!(values(bqm_copy), Vector{Int}[], var)
    bqm_copy
end

+(var::Real, bqm::BinaryOptimizationModel) = +(bqm, var)
-(bqm::BinaryOptimizationModel, var::Real) = +(bqm, -var)

## sum of bqm
function _bqm_add(bqm1::S, bqm2::S) where S<:BinaryOptimizationModel{<:Real}
    bqm_result = S(vars(bqm1), copy(values(bqm1)))
    for (key, val) = values(bqm2)
        _dict_add!(values(bqm_result), key, val)
    end
    bqm_result
end

function +(bqm1::BinaryOptimizationModel, bqm2::BinaryOptimizationModel)
    if length(values(bqm1)) > length(values(bqm2))
        _bqm_add(bqm1, bqm2)
    else
        _bqm_add(bqm2, bqm1)
    end
end

-(bqm1::BinaryOptimizationModel, bqm2::BinaryOptimizationModel) = bqm1 + (-bqm2)
-(var::Float64, bqm::BinaryOptimizationModel) = var + (-bqm)

## Multiplication of hobos
function *(bqm1::S, bqm2::S) where S<:BinaryOptimizationModel{T} where T<:Real
    @assert vars(bqm1) == vars(bqm2)
    bqm_result = S(vars(bqm1))
    for (key1, val1) = values(bqm1), (key2, val2) = values(bqm2)
        _dict_add!(values(bqm_result), clean_key(S, vcat(key1, key2)), val1*val2)
    end
    bqm_result
end

##
function order(bqm::BinaryOptimizationModel)
    if !isempty(values(bqm))
        maximum(length.(keys(values(bqm))))
    else
        0
    end
end

isquadratic(bqm::BinaryOptimizationModel) = order(bqm) <= 2


text:
module STMO

using Plots, LinearAlgebra

include("helpers.jl")
export myred, mygreen, myyellow, myblue, myblack, myorange, mycolors
export trace, notrack, PathTrack, nsteps, Tracker
export path, path!, plotobj, plotobj!, colorscatter
export dist, hamming
export EdgeList, WeightedEdgeList, Vertices, AdjList, edges2adjlist, adjlist2edges
export nvertices, isconnected, vertices

include("quadratic.jl")
export fquad, quadratic
export Quadratic, SignalRecovery
#export solve_quadratic, quadratic_ls, gradient_descent, plot_quadratic
#export generate_noisy_measurements, make_bookkeeping, make_connection_matrix, signalfun

include("unconstrained.jl")
export Unconstrained
#export gradient_descent, coordinate_descent, newtons_method, backtracking_line_search

include("optimaltransport.jl")
export OptimalTransport

include("mst.jl")
export MST

include("shortestpath.jl")
export ShortestPath

include("testfuns.jl")
export TestFuns

include("tickettoride.jl")
export TicketToRide

include("words.jl")
export Words

include("cancer.jl")
export CancerData

include("tsp.jl")
export TSP

end # module


text:
using MSM
using Documenter

DocMeta.setdocmeta!(MSM, :DocTestSetup, :(using MSM); recursive=true)

makedocs(;
    modules=[MSM],
    authors="banachtech <balaji@banach.tech> and contributors",
    repo="https://github.com/banachtech/MSM.jl/blob/{commit}{path}#{line}",
    sitename="MSM.jl",
    format=Documenter.HTML(;
        prettyurls=get(ENV, "CI", "false") == "true",
        canonical="https://banachtech.github.io/MSM.jl",
        assets=String[],
    ),
    pages=[
        "Home" => "index.md",
    ],
)

deploydocs(;
    repo="github.com/banachtech/MSM.jl",
)


text:
abstract type AbstractField end

struct ScalarField{D,G,U} <: AbstractField
    data::D
    grid::G
    utils::U
end

function ScalarField(data, grid)
    ℱ = plan_fft!(data)
    ℱ⁻¹ = plan_ifft!(data)
    utils = (; forward=ℱ, backward=ℱ⁻¹)
    return ScalarField(data, grid, utils)
end

struct VectorField{D, G, U} <: AbstractField
    data::D
    grid::G
    utils::U
end

struct StateField{D, G, U} <: AbstractField
    data::D
    grid::G
    utils::U
end

max_stars_repo_name: sandreza/FourierCore

max_issues_repo_name: sandreza/FourierCore

max_forks_repo_name: sandreza/FourierCore


text:
# # 2D process with delayed-alpha mutual interactions and exp self interaction

#=
Here I consider two neurons with a delayed interaction,
that self-interact without delay.
=#

# ## Initialization

using ProgressMeter
using LinearAlgebra,Statistics,StatsBase,Distributions
using Plots,NamedColors ; theme(:default)
using SparseArrays
using FFTW
using Random
Random.seed!(0)

using HawkesSimulator; global const H = HawkesSimulator

##

# ## Define the interaction kernel
# The kernel is a delayed-alpha, see that example file for definition and
# details.
#
# Important : the autaptic weights (diagonal of `wmat`) are set to zero
# if not, there would also be a self-delayed-interaction.

mytau = 0.6
mydelay = 1.0
mywmat = [ 0.0   0.1
           1.8   0.0 ]

myinput = [0.2,0.1]

pop = H.PopulationAlphaDelay(mytau,mydelay)
pops = H.PopulationState(pop,myinput)

# ## Define the autapses
# the autaptic weights are `mywautaps`
# the self-interaction kernel is `popautaps`

mywautaps = [0.5,0.2]
popautaps = H.PopulationExp(0.4)
autapses_pop = H.Autapses(mywautaps,popautaps)

# ## Build the network and run simulation

netw = H.InputNetwork(pops,[pops,],[mywmat,],autapses_pop)

function run_simulation!(netw,nspikes)
  t_now = 0.0
  H.reset!(netw) # clear spike trains etc
  @showprogress 1.0 "Running Hawkes process..." for k in 1:nspikes
    t_now = H.dynamics_step!(t_now,[netw,])
    if k%1_000 == 0
      H.clear_trains!(netw.postpops,50.0) # longer time horizon, because of delay
    end
  end
  H.clear_trains!(netw.postpops,-1.0);
  return t_now
end

n_spikes = 100_000
Tmax = run_simulation!(ntw,n_spikes);

@info "Total duration $(round(Tmax;digits=1))"

# ## Compare rate , simulation vs analytic

ratenum = H.numerical_rates(pops)

# analytic rate from fourier Eq between 6 and  7 in Hawkes
# note that I need to add autapses to `G0` separately

ratefou = let G0 =  mywmat .* H.interaction_kernel_fourier(0,pop)
  G0 .+=  diagm( 0 => mywautaps .* H.interaction_kernel_fourier(0,popautaps) )
  inv(I-G0)*myinput|> real
end

@info "Mean rate -  simulation $(round.(ratenum;digits=2)), analytic  $(round.(ratefou;digits=2))"

##
# ## Covariance density
# As in previous examples, compute from simulaton and then analytic

myspikes_all = pops.trains_history

mydt = 0.1
myτmax = 30.0
mytaus = H.get_times(mydt,myτmax)
ntaus = length(mytaus)
cov_num = H.covariance_density_numerical(myspikes_all,mydt,myτmax)

# plot all numerical covariance densities

function doplot()
  plt = plot(xlabel="time delay (s)",ylabel="Covariance density",title="simulation")
  for i in 1:2, j in 1:2
    plot!(plt,mytaus[2:end-1],cov_num[i,j,2:end-1] ; linewidth = 2,label="cov $i - $j")
  end
  return plt
end
doplot()

# the bumps reflect the delayed interactions

function four_high_res(dt::Real,Tmax::Real)
  k1 = 2
  k2 = 0.01
  myτmax,mydt = Tmax * k1, dt*k2
  mytaus = H.get_times(mydt,myτmax)
  nkeep = div(length(mytaus),k1)
  myfreq = H.get_frequencies_centerzero(mydt,myτmax)
  G_omega = map(mywmat) do w
    ifftshift( w .* H.interaction_kernel_fourier.(myfreq,Ref(pop)))
  end
  for k in 1:size(mywmat,1) # add autapses
    waut = mywautaps[k]
    Gaut = ifftshift(waut .*  H.interaction_kernel_fourier.(myfreq,Ref(popautaps)))
    G_omega[k,k] .+= Gaut
  end
  D = Diagonal(ratefou)
  M = Array{ComplexF64}(undef,2,2,length(myfreq))
  Mt = similar(M,Float64)
  for i in eachindex(myfreq)
    G = getindex.(G_omega,i)
    M[:,:,i] = (I-G)\D/(I-G')
  end
  for i in 1:2,j in 1:2
    Mt[i,j,:] = real.(ifft(M[i,j,:]))
    Mt[i,j,2:end] ./= mydt
  end
  return mytaus[1:nkeep],Mt[:,:,1:nkeep]
end

(taush,Cfou)=four_high_res(mydt,myτmax)

function oneplot(i,j)
  plt = plot(xlabel="time delay (s)",ylabel="Covariance density",title="cov $i - $j")
  plot!(plt,mytaus[2:end],cov_num[i,j,2:end] ; linewidth = 3)
  return plot!(plt,taush[2:end],Cfou[i,j,2:end]; linestyle=:dash, linewidth=3)
end

# 1
oneplot(1,1)
# 2
oneplot(1,2)
# 3
oneplot(2,1)
# 4
oneplot(2,2)

# Perfect match!

##
# ## Visualize raster plot
# the raster plot shows some correlation between the neural activities
# neuron one (lower row) excites neuron two (upper row)

function doplot(tlims = (3000.,3080.) )
  plt=plot()
  for (is,_train) in enumerate(myspikes_all)
    train = filter(t-> tlims[1]< t < tlims[2],_train)
    nspk = length(train)
    scatter!(plt,train,fill(0.1*(is),nspk),markersize=30,
      markercolor=:black,markershape=:vline,leg=false)
  end
  return plot!(plt,ylims=(0.0,0.3))
end

doplot()

# **THE END**

# Literate.markdown("examples/2d_delay_autapses.jl","docs/src";documenter=true,repo_root_url="https://github.com/dylanfesta/HawkesSimulator.jl/blob/master") #src


text:
function get_individual(survey, outfile = "temp.xlsx")
	url = url_individual(survey)
	urldownload(url, parser = identity, save_raw = outfile)

	df = DataFrame(XLSX.readtable(outfile, survey)...)
	df = ifelse.(df .== "#N/A", missing, df)
	df.ID = convert.(Int, df.ID)

	outfile == "temp.xlsx" ? rm(outfile) : nothing
	return df
end


text:
#---------------
# copying matrix
#---------------
function mycopy_mat(MATin,MATout,nz,nr)
    @inbounds Threads.@threads for j=1:nr
        for k=1:nz
            MATout[k,j]=MATin[k,j]
        end
    end
end

#------------------------------
# finding acoustic or elastic media,
# and assigning flags where
# vf should be zero
#------------------------------
function get_Flag_vf_zero(Flag_AC,Flag_E,nr,nz)
    #when FlagAC==1 or Flag_E==1, then set Flag_vf_zero==1
    Flag_vf_zero=zeros(nz,nr)
    @inbounds Threads.@threads for j=1:nr
        for k=1:nz
            if(Flag_AC[k,j]==1 || Flag_E[k,j]==1)
                Flag_vf_zero[k,j]=1
            end
        end
    end
    return Flag_vf_zero
end

#------------------------------
# Checking stability of FD modeling
#------------------------------
function check_stability01(dt,dr,Vp,m)
#From Randall 1991, cdt/dr<=1/sqrt(2)*(1+m^2dr^2/(8*rj^2))^(-1/2)
#dr=dz, rj=j*dr
    #Rough version : V*dt/dr<=1/sqrt(2), V-->maximum Vp
    Vp_max=maximum(Vp)
    if (Vp_max*dt/dr <= 1/sqrt(2))
        println("Stability OK.")
    else
        error("########Stability NG.######")
    end

end

#------------------------------------------------------------
# Following are not used, but kept for future extention
#------------------------------------------------------------
#--calculate divergence of pariticle velocity (assuming axisymetric)
function wave_div01!(divV,nr,nz,dr,dz,vr,vz)
   #div(u)=1/r*dur/dr+duz/dz: location of vphi,tii (Mittet)
#@inbounds Threads.@threads
   for j=3:nr-2
      for k=3:nz-2
         r_now=(j-1)*dr #position of tii (nonaxis to avoid singularity)
         vr_f1,vr_b1=vr[k,j],vr[k,j-1]
         vr_f2,vr_b2=vr[k,j+1],vr[k,j-2]
         vz_f1,vz_b1=vz[k,j],vz[k-1,j]
         vz_f2,vz_b2=vz[k+1,j],vz[k-2,j]
         drvr=9.0/8.0*(vr_f1-vr_b1)-1.0/24.0*(vr_f2-vr_b2)
         drvr=drvr/dr
         dzvz=9.0/8.0*(vz_f1-vz_b1)-1.0/24.0*(vz_f2-vz_b2)
         dzvz=dzvz/dz
         divV[k,j]=1.0/r_now*drvr+dzvz
      end
   end
end

#--calculate curl of pariticle velocity (assuming axisymetric)
function wave_curl01!(curV,nr,nz,dr,dz,vr,vz)
   #curl(u)=dur/dz-duz/dr: location of trz (Mittet)
#@inbounds Threads.@threads
   for j=3:nr-2
      for k=3:nz-2
         r_now=(j-1)*dr+dr/2.0 #position of trz
         vr_f1,vr_b1=vr[k+1,j],vr[k,j]
         vr_f2,vr_b2=vr[k+2,j],vr[k-1,j]
         vz_f1,vz_b1=vz[k,j+1],vz[k,j]
         vz_f2,vz_b2=vz[k,j+2],vz[k,j-1]
         dzvr=9.0/8.0*(vr_f1-vr_b1)-1.0/24.0*(vr_f2-vr_b2)
         dzvr=dzvr/dz
         drvz=9.0/8.0*(vz_f1-vz_b1)-1.0/24.0*(vz_f2-vz_b2)
         drvz=drvz/dr
         curV[k,j]=dzvr-drvz
      end
   end
end

function init_snap_sparse(nt,nz,nr)
   #---Initializing snapshots (sparse version)
   nskip=100
   nr_sp=0
   for ir=1:10:nr
      nr_sp=nr_sp+1
   end
   nz_sp=0
   for iz=1:10:nz
      nz_sp=nz_sp+1
   end

   snapshots_vr,snapshots_vz,nsnap,itvec_snap=init_snapshots_v(nskip,nt,nz_sp,nr_sp)
   snapshots_trr,nsnap,itvec_snap=init_snapshots_t(nskip,nt,nz_sp,nr_sp)
   return nskip,snapshots_trr,snapshots_vr,snapshots_vz,nsnap,itvec_snap
end


text:
function set_lane_changing_pts!(model::MultiPtsDriver, scene::Union{Scene,Frame{Entity{VehicleState, BicycleModel, Int}}}, roadway::Roadway, ego_index::Int)
    # direction is like -1: turn to right lane, 1: turn to left lane
    #println("enter function")

    veh = scene[ego_index]
    x = veh.state.posG.x
    y = veh.state.posG.y
    θ = veh.state.posG.θ
    steer=model.siji.steer
    steerdotmax=model.steerdotmax
    steermax=model.steermax
    accmax=model.accmax
    v = veh.state.v
    roadind=veh.state.posF.roadind
    Δt = model.Δt

    #determine transitionStep
    #if it is current lane
    #first transition
    #second transition
    #smoothing
    laneNum=roadind.tag.lane
    transitionStep=0
    v_desire=model.v_desire

    pos=veh.state.posG
    transitionMagnitude=0.2
    if model.laneNum_desire>laneNum #turn to left lane
        transitionStep=transitionMagnitude
        roadind1=roadind
        frenet=Frenet(roadind1, roadway[roadind1].s,veh.state.posF.t,veh.state.posF.ϕ)
        i=1
        while frenet.t+transitionStep<DEFAULT_LANE_WIDTH/2 && i<=20

            v_desire_n=v_desire+model.acc_desire*(Δt*2)
            v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
            Δs=(v_desire+v_desire_n)*(Δt*2)/2
            v_desire=v_desire_n

            if i==6
                model.v_desire=v_desire
            end

            roadind1=move_along(roadind1, roadway, Δs)
            frenet=Frenet(roadind1, roadway[roadind1].s,frenet.t+transitionStep,veh.state.posF.ϕ)
            pos=get_posG(frenet,roadway)
            model.Pts[1,i]=pos.x
            model.Pts[2,i]=pos.y
            i=i+1
        end
        if frenet.t+transitionStep>=DEFAULT_LANE_WIDTH/2 && i<=20

            v_desire_n=v_desire+model.acc_desire*(Δt*2)
            v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
            Δs=(v_desire+v_desire_n)*(Δt*2)/2
            v_desire=v_desire_n

            if i==6
                model.v_desire=v_desire
            end

            lane = roadway[roadind1.tag]
            lane_desire = roadway[LaneTag(lane.tag.segment, lane.tag.lane + 1)]
            posG=pos
            roadproj = proj(posG, lane_desire, roadway)
            roadind2=RoadIndex(roadproj)
            roadind2=move_along(roadind2, roadway, Δs)
            frenet=Frenet(roadind2, roadway[roadind2].s,-DEFAULT_LANE_WIDTH/2+frenet.t+transitionStep-DEFAULT_LANE_WIDTH/2,veh.state.posF.ϕ)
            pos=get_posG(frenet,roadway)
            model.Pts[1,i]=pos.x
            model.Pts[2,i]=pos.y
            i=i+1
            while frenet.t+transitionStep<=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind2=move_along(roadind2, roadway, Δs)
                frenet=Frenet(roadind2, roadway[roadind2].s,frenet.t+transitionStep,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
            while frenet.t+transitionStep>=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind2=move_along(roadind2, roadway, Δs)
                frenet=Frenet(roadind2, roadway[roadind2].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        end

    elseif model.laneNum_desire<laneNum #turn to right lane
        transitionStep=-transitionMagnitude
        roadind1=roadind
        frenet=Frenet(roadind1, roadway[roadind1].s,veh.state.posF.t,veh.state.posF.ϕ)
        i=1
        while frenet.t+transitionStep>-DEFAULT_LANE_WIDTH/2 && i<=20

            v_desire_n=v_desire+model.acc_desire*(Δt*2)
            v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
            Δs=(v_desire+v_desire_n)*(Δt*2)/2
            v_desire=v_desire_n

            if i==6
                model.v_desire=v_desire
            end

            roadind1=move_along(roadind1, roadway, Δs)
            frenet=Frenet(roadind1, roadway[roadind1].s,frenet.t+transitionStep,veh.state.posF.ϕ)
            pos=get_posG(frenet,roadway)
            model.Pts[1,i]=pos.x
            model.Pts[2,i]=pos.y
            i=i+1
        end
        if frenet.t+transitionStep<=DEFAULT_LANE_WIDTH/2 && i<=20

            v_desire_n=v_desire+model.acc_desire*(Δt*2)
            v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
            Δs=(v_desire+v_desire_n)*(Δt*2)/2
            v_desire=v_desire_n

            if i==6
                model.v_desire=v_desire
            end

            lane = roadway[roadind1.tag]
            lane_desire = roadway[LaneTag(lane.tag.segment, lane.tag.lane - 1)]
            posG=pos
            roadproj = proj(posG, lane_desire, roadway)
            roadind2=RoadIndex(roadproj)
            roadind2=move_along(roadind2, roadway, Δs)
            frenet=Frenet(roadind2, roadway[roadind2].s,DEFAULT_LANE_WIDTH/2+(frenet.t+transitionStep+DEFAULT_LANE_WIDTH/2),veh.state.posF.ϕ)
            pos=get_posG(frenet,roadway)
            model.Pts[1,i]=pos.x
            model.Pts[2,i]=pos.y
            i=i+1
            while frenet.t+transitionStep>=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind2=move_along(roadind2, roadway, Δs)
                frenet=Frenet(roadind2, roadway[roadind2].s,frenet.t+transitionStep,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
            while frenet.t+transitionStep<=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind2=move_along(roadind2, roadway, Δs)
                frenet=Frenet(roadind2, roadway[roadind2].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        end

    else #on the same lane try keep driving on the center lane
        if veh.state.posF.t>0
            transitionStep=-transitionMagnitude
            roadind1=roadind
            frenet=Frenet(roadind1, roadway[roadind1].s,veh.state.posF.t,veh.state.posF.ϕ)
            i=1
            while frenet.t+transitionStep>0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along(roadind1, roadway, Δs)
                frenet=Frenet(roadind1, roadway[roadind1].s,frenet.t+transitionStep,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
            while frenet.t+transitionStep<=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along(roadind1, roadway, Δs)
                frenet=Frenet(roadind1, roadway[roadind1].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        elseif veh.state.posF.t<0
            transitionStep=transitionMagnitude
            roadind1=roadind
            frenet=Frenet(roadind1, roadway[roadind1].s,veh.state.posF.t,veh.state.posF.ϕ)
            i=1
            while frenet.t+transitionStep<0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along(roadind1, roadway, Δs)
                frenet=Frenet(roadind1, roadway[roadind1].s,frenet.t+transitionStep,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
            while frenet.t+transitionStep>=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along(roadind1, roadway, Δs)
                frenet=Frenet(roadind1, roadway[roadind1].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        elseif veh.state.posF.t==0
            roadind1=roadind
            frenet=Frenet(roadind1, roadway[roadind1].s,0.0,veh.state.posF.ϕ)
            i=1
            while i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along(roadind1, roadway, Δs)
                frenet=Frenet(roadind1, roadway[roadind1].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        end
    end
    model.index=1
    model.subindex=0
    #getAccSteer!(model,scene,roadway,ego_index)
    model
end

#########

function set_lane_changing_pts_with_direction!(model::MultiPtsDriver, scene::Union{Scene,Frame{Entity{VehicleState, BicycleModel, Int}}}, roadway::Roadway, ego_index::Int, turning_direction::Int = 1)

    veh = scene[ego_index]
    x = veh.state.posG.x
    y = veh.state.posG.y
    θ = veh.state.posG.θ
    steer=model.siji.steer
    steerdotmax=model.steerdotmax
    steermax=model.steermax
    accmax=model.accmax
    v = veh.state.v
    roadind=veh.state.posF.roadind
    Δt = model.Δt

    #determine transitionStep
    #if it is current lane
    #first transition
    #second transition
    #smoothing
    laneNum=roadind.tag.lane
    transitionStep=0
    v_desire=model.v_desire

    pos=veh.state.posG
    transitionMagnitude=0.2
    if model.laneNum_desire>laneNum #turn to left lane
        transitionStep=transitionMagnitude
        roadind1=roadind
        frenet=Frenet(roadind1, roadway[roadind1].s,veh.state.posF.t,veh.state.posF.ϕ)
        i=1
        while frenet.t+transitionStep<DEFAULT_LANE_WIDTH/2 && i<=20

            v_desire_n=v_desire+model.acc_desire*(Δt*2)
            v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
            Δs=(v_desire+v_desire_n)*(Δt*2)/2
            v_desire=v_desire_n

            if i==6
                model.v_desire=v_desire
            end

            roadind1=move_along_with_direction(roadind1, roadway, Δs,direction=turning_direction)
            frenet=Frenet(roadind1, roadway[roadind1].s,frenet.t+transitionStep,veh.state.posF.ϕ)
            pos=get_posG(frenet,roadway)
            model.Pts[1,i]=pos.x
            model.Pts[2,i]=pos.y
            i=i+1
        end
        if frenet.t+transitionStep>=DEFAULT_LANE_WIDTH/2 && i<=20

            v_desire_n=v_desire+model.acc_desire*(Δt*2)
            v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
            Δs=(v_desire+v_desire_n)*(Δt*2)/2
            v_desire=v_desire_n

            if i==6
                model.v_desire=v_desire
            end

            lane = roadway[roadind1.tag]
            lane_desire = roadway[LaneTag(lane.tag.segment, lane.tag.lane + 1)]
            posG=pos
            roadproj = proj(posG, lane_desire, roadway)
            roadind2=RoadIndex(roadproj)
            roadind2=move_along_with_direction(roadind2, roadway, Δs,direction=turning_direction)
            frenet=Frenet(roadind2, roadway[roadind2].s,-DEFAULT_LANE_WIDTH/2+frenet.t+transitionStep-DEFAULT_LANE_WIDTH/2,veh.state.posF.ϕ)
            pos=get_posG(frenet,roadway)
            model.Pts[1,i]=pos.x
            model.Pts[2,i]=pos.y
            i=i+1
            while frenet.t+transitionStep<=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind2=move_along_with_direction(roadind2, roadway, Δs,direction=turning_direction)
                frenet=Frenet(roadind2, roadway[roadind2].s,frenet.t+transitionStep,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
            while frenet.t+transitionStep>=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind2=move_along_with_direction(roadind2, roadway, Δs,direction=turning_direction)
                frenet=Frenet(roadind2, roadway[roadind2].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        end

    elseif model.laneNum_desire<laneNum #turn to right lane
        transitionStep=-transitionMagnitude
        roadind1=roadind
        frenet=Frenet(roadind1, roadway[roadind1].s,veh.state.posF.t,veh.state.posF.ϕ)
        i=1
        while frenet.t+transitionStep>-DEFAULT_LANE_WIDTH/2 && i<=20

            v_desire_n=v_desire+model.acc_desire*(Δt*2)
            v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
            Δs=(v_desire+v_desire_n)*(Δt*2)/2
            v_desire=v_desire_n

            if i==6
                model.v_desire=v_desire
            end

            roadind1=move_along_with_direction(roadind1, roadway, Δs,direction=turning_direction)
            frenet=Frenet(roadind1, roadway[roadind1].s,frenet.t+transitionStep,veh.state.posF.ϕ)
            pos=get_posG(frenet,roadway)
            model.Pts[1,i]=pos.x
            model.Pts[2,i]=pos.y
            i=i+1
        end
        if frenet.t+transitionStep<=DEFAULT_LANE_WIDTH/2 && i<=20

            v_desire_n=v_desire+model.acc_desire*(Δt*2)
            v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
            Δs=(v_desire+v_desire_n)*(Δt*2)/2
            v_desire=v_desire_n

            if i==6
                model.v_desire=v_desire
            end

            lane = roadway[roadind1.tag]
            lane_desire = roadway[LaneTag(lane.tag.segment, lane.tag.lane - 1)]
            posG=pos
            roadproj = proj(posG, lane_desire, roadway)
            roadind2=RoadIndex(roadproj)
            roadind2=move_along_with_direction(roadind2, roadway, Δs,direction=turning_direction)
            frenet=Frenet(roadind2, roadway[roadind2].s,DEFAULT_LANE_WIDTH/2+(frenet.t+transitionStep+DEFAULT_LANE_WIDTH/2),veh.state.posF.ϕ)
            pos=get_posG(frenet,roadway)
            model.Pts[1,i]=pos.x
            model.Pts[2,i]=pos.y
            i=i+1
            while frenet.t+transitionStep>=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind2=move_along_with_direction(roadind2, roadway, Δs,direction=turning_direction)
                frenet=Frenet(roadind2, roadway[roadind2].s,frenet.t+transitionStep,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
            while frenet.t+transitionStep<=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind2=move_along_with_direction(roadind2, roadway, Δs,direction=turning_direction)
                frenet=Frenet(roadind2, roadway[roadind2].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        end

    else #on the same lane try keep driving on the center lane
        if veh.state.posF.t>0
            transitionStep=-transitionMagnitude
            roadind1=roadind
            frenet=Frenet(roadind1, roadway[roadind1].s,veh.state.posF.t,veh.state.posF.ϕ)
            i=1
            while frenet.t+transitionStep>0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along_with_direction(roadind1, roadway, Δs,direction=turning_direction)
                frenet=Frenet(roadind1, roadway[roadind1].s,frenet.t+transitionStep,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
            while frenet.t+transitionStep<=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along_with_direction(roadind1, roadway, Δs,direction=turning_direction)
                frenet=Frenet(roadind1, roadway[roadind1].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        elseif veh.state.posF.t<0
            transitionStep=transitionMagnitude
            roadind1=roadind
            frenet=Frenet(roadind1, roadway[roadind1].s,veh.state.posF.t,veh.state.posF.ϕ)
            i=1
            while frenet.t+transitionStep<0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along_with_direction(roadind1, roadway, Δs,direction=turning_direction)
                frenet=Frenet(roadind1, roadway[roadind1].s,frenet.t+transitionStep,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
            while frenet.t+transitionStep>=0 && i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along_with_direction(roadind1, roadway, Δs,direction=turning_direction)
                frenet=Frenet(roadind1, roadway[roadind1].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        elseif veh.state.posF.t==0
            roadind1=roadind
            frenet=Frenet(roadind1, roadway[roadind1].s,0.0,veh.state.posF.ϕ)
            i=1
            while i<=20

                v_desire_n=v_desire+model.acc_desire*(Δt*2)
                v_desire_n=clamp(v_desire_n,model.v_min,model.v_max)
                Δs=(v_desire+v_desire_n)*(Δt*2)/2
                v_desire=v_desire_n

                if i==6
                    model.v_desire=v_desire
                end

                roadind1=move_along_with_direction(roadind1, roadway, Δs,direction=turning_direction)
                frenet=Frenet(roadind1, roadway[roadind1].s,0.0,veh.state.posF.ϕ)
                pos=get_posG(frenet,roadway)
                model.Pts[1,i]=pos.x
                model.Pts[2,i]=pos.y
                i=i+1
            end
        end
    end
    model.index=1
    model.subindex=0
    #getAccSteer!(model,scene,roadway,ego_index)
    model
end

max_stars_repo_name: sisl/UrbanDrivingSimulation

max_issues_repo_name: sisl/UrbanDrivingSimulation

max_forks_repo_name: sisl/UrbanDrivingSimulation


text:
# Note that this script can accept some limited command-line arguments, run
# `julia build_tarballs.jl --help` to see a usage message.
using BinaryBuilder, Pkg

name = "Mineos"
version = v"1.0.1" # Artificial version number for bumping Julia compat (and expanding platforms)

# Collection of sources required to build Mineos
sources = [
    GitSource("https://github.com/anowacki/mineos.git", "e2558b486d7656ef112608a8776643da66dc87cf"),
]

# Bash recipe for building across all platforms
script = raw"""
cd ${WORKSPACE}/srcdir/mineos

# Fix clang error 'error: non-void function * should return a value [-Wreturn-type]'
if [[ "${target}" == *-freebsd* ]] || [[ "${target}" == *-apple-* ]]; then
    export CFLAGS="-Wno-return-type"
fi

# Fix issue due to GCC 10+. Only necessary on aarch64-apple as others use earlier GCC
#     Error: Rank mismatch between actual argument at (1) and actual argument at (2) (scalar and rank-1)
if [[ "${target}" == aarch64-apple-* ]]; then
    export FFLAGS="-std=legacy"
fi

./configure --prefix=${prefix} --build=${MACHTYPE} --host=${target} --disable-doc
make

make install
"""

# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
platforms = expand_gfortran_versions(supported_platforms())

# The products that we will ensure are always built
products = [
    ExecutableProduct("simpledit", :simpledit),
    ExecutableProduct("endi", :endi),
    ExecutableProduct("eigcon", :eigcon),
    ExecutableProduct("eigen2asc", :eigen2asc),
    ExecutableProduct("green", :green),
    ExecutableProduct("syndat", :syndat),
    ExecutableProduct("minos_bran", :minos_bran),
    ExecutableProduct("cucss2sac", :cucss2sac)
]

# Dependencies that must be installed before this package can be built
dependencies = [
    Dependency(PackageSpec(name="CompilerSupportLibraries_jll", uuid="e66e0078-7015-5450-92f7-15fbd957f2ae"))
]

# Build the tarballs, and possibly a `build.jl` as well.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies; julia_compat = "1.6")

max_stars_repo_name: jakebolewski/Yggdrasil

max_issues_repo_name: jakebolewski/Yggdrasil

max_forks_repo_name: jakebolewski/Yggdrasil