text:
__precompile__()

module Compat

using Base.Meta

"""Get just the function part of a function declaration."""
withincurly(ex) = isexpr(ex, :curly) ? ex.args[1] : ex

if VERSION < v"0.6.0-dev.2043"
    Base.take!(t::Task) = consume(t)
end

function rewrite_show(ex)
    if isexpr(ex, :call)
        Expr(:call, rewrite_show(ex.args[1]), ex.args[2:end]...)
    elseif isexpr(ex, :curly)
        Expr(:curly, rewrite_show(ex.args[1]), ex.args[2:end]...)
    else
        :(Base.writemime)
    end
end

function rewrite_dict(ex)
    length(ex.args) == 1 && return ex

    f = ex.args[1]
    if isexpr(f, :curly)
        newex = Expr(:typed_dict, :($(f.args[2])=>$(f.args[3])))
    else
        newex = Expr(:dict)
    end

    for i = 2:length(ex.args)
        pair = ex.args[i]
        !isexpr(pair, :(=>)) && return ex
        push!(newex.args, pair)
    end
    newex
end

function rewrite_ordereddict(ex)
    length(ex.args) == 1 && return ex

    f = ex.args[1]
    newex = Expr(:call, f, :[])

    for i = 2:length(ex.args)
        pair = ex.args[i]
        !isexpr(pair, :(=>)) && return ex
        push!(newex.args[2].args, Expr(:tuple, pair.args...))
    end

    newex
end

# rewrites all subexpressions of the form `a => b` to `(a, b)`
function rewrite_pairs_to_tuples!(expr::Expr)
    if expr.head == :(=>)
        expr.head = :tuple
    end
    for subexpr in expr.args
        isa(subexpr, Expr) && rewrite_pairs_to_tuples!(subexpr)
    end
    return expr
end

function is_quote_symbol(ex::ANY, val::Symbol)
    if isa(ex, QuoteNode)
        return (ex::QuoteNode).value === val
    elseif isa(ex, Expr)
        ex = ex::Expr
        return ex.head === :quote && length(ex.args) == 1 && ex.args[1] === val
    end
    return false
end

is_index_style(ex::Expr) = ex == :(Compat.IndexStyle) || ex == :(Base.IndexStyle) ||
    (ex.head == :(.) && (ex.args[1] == :Compat || ex.args[1] == :Base) &&
         ex.args[2] == Expr(:quote, :IndexStyle))

is_index_style(arg) = false

# rewrites accesses to IOContext dicts
function rewrite_iocontext!(expr::Expr)
    args = expr.args
    nargs = length(args)
    if nargs == 4 && expr.head === :call && args[1] === :get && args[4] === false
        key = args[3]
        if is_quote_symbol(key, :limit) || is_quote_symbol(key, :compact)
            if VERSION >= v"0.5.0-dev+1936" && VERSION < v"0.5.0-dev+4305"
                args[1] = :(Base.limit_output)
                deleteat!(args, 3:4)
            elseif VERSION < v"0.5.0-dev+1936"
                expr.head = :quote
                args[1] = false
                deleteat!(args, 3:4)
            end
        elseif is_quote_symbol(key, :multiline)
            if VERSION < v"0.5.0-dev+4305"
                expr.head = :quote
                args[1] = false
                deleteat!(args, 3:4)
            end
        end
    end
end

# JuliaLang/julia#10543
if !isdefined(Base, :tryparse)
    function tryparse{T}(::Type{T}, args...)
        try
            Nullable(Base.parse(T, args...))
        catch
            Nullable{T}()
        end
    end
end

import Base.unsafe_convert

function new_style_call_overload(ex::Expr)
    Base.depwarn("new_style_call_overload is deprecated.", :new_style_call_overload)
    false
end

istopsymbol(ex, mod, sym) = ex in (sym, Expr(:(.), mod, Expr(:quote, sym)))

if VERSION < v"0.6.0-dev.2782"
    function new_style_typealias(ex::ANY)
        isexpr(ex, :(=)) || return false
        ex = ex::Expr
        return length(ex.args) == 2 && isexpr(ex.args[1], :curly)
    end
else
    new_style_typealias(ex) = false
end

function _compat(ex::Expr)
    if ex.head === :call
        f = ex.args[1]
        if VERSION < v"0.6.0-dev.826" && length(ex.args) == 3 && # julia#18510
                istopsymbol(withincurly(ex.args[1]), :Base, :Nullable)
            ex = Expr(:call, f, ex.args[2], Expr(:call, :(Compat._Nullable_field2), ex.args[3]))
        end
    elseif ex.head === :curly
        f = ex.args[1]
        if ex == :(Ptr{Void})
            # Do not change Ptr{Void} to Ptr{Nothing}: 0.4.0-dev+768
            return ex
        elseif VERSION < v"0.6.0-dev.2575" #20414
            ex = Expr(:curly, map(a -> isexpr(a, :call, 2) && a.args[1] == :(<:) ?
                                  :($TypeVar($(QuoteNode(gensym(:T))), $(a.args[2]), false)) :
                                  isexpr(a, :call, 2) && a.args[1] == :(>:) ?
                                  :($TypeVar($(QuoteNode(gensym(:T))), $(a.args[2]), $Any, false)) : a,
                                  ex.args)...)
        end
    elseif ex.head === :macrocall
        f = ex.args[1]
    elseif ex.head === :quote && isa(ex.args[1], Symbol)
        # Passthrough
        return ex
    elseif new_style_typealias(ex)
        ex.head = :typealias
    elseif ex.head === :const && length(ex.args) == 1 && new_style_typealias(ex.args[1])
        ex = ex.args[1]::Expr
        ex.head = :typealias
    end
    if VERSION < v"0.6.0-dev.2840"
        if ex.head == :(=) && isa(ex.args[1], Expr) && ex.args[1].head == :call
            a = ex.args[1].args[1]
            if is_index_style(a)
                ex.args[1].args[1] = :(Base.linearindexing)
            elseif isa(a, Expr) && a.head == :curly
                if is_index_style(a.args[1])
                    ex.args[1].args[1].args[1] = :(Base.linearindexing)
                end
            end
        end
    end
    return Expr(ex.head, map(_compat, ex.args)...)
end

_compat(ex) = ex

function _get_typebody(ex::Expr)
    args = ex.args
    if ex.head !== :type || length(args) != 3 || args[1] !== true
        throw(ArgumentError("Invalid usage of @compat: $ex"))
    end
    name = args[2]
    if !isexpr(args[3], :block)
        throw(ArgumentError("Invalid type declaration: $ex"))
    end
    body = (args[3]::Expr).args
    filter!(body) do e
        if isa(e, LineNumberNode) || isexpr(e, :line)
            return false
        end
        return true
    end
    return name, body
end

function _compat_primitive(typedecl)
    name, body = _get_typebody(typedecl)
    if length(body) != 1
        throw(ArgumentError("Invalid primitive type declaration: $typedecl"))
    end
    return Expr(:bitstype, body[1], name)
end

function _compat_abstract(typedecl)
    name, body = _get_typebody(typedecl)
    if length(body) != 0
        throw(ArgumentError("Invalid abstract type declaration: $typedecl"))
    end
    return Expr(:abstract, name)
end

macro compat(ex...)
    if VERSION < v"0.6.0-dev.2746" && length(ex) == 2 && ex[1] === :primitive
        return esc(_compat_primitive(ex[2]))
    elseif length(ex) != 1
        throw(ArgumentError("@compat called with wrong number of arguments: $ex"))
    elseif (VERSION < v"0.6.0-dev.2746" && isexpr(ex[1], :abstract) &&
            length(ex[1].args) == 1 && isexpr(ex[1].args[1], :type))
        # This can in principle be handled in nested case but we do not
        # do that to be consistent with primitive types.
        return esc(_compat_abstract(ex[1].args[1]))
    end
    esc(_compat(ex[1]))
end

export @compat, @inline, @noinline

import Base.@irrational

import Base: remotecall, remotecall_fetch, remotecall_wait, remote_do

import Base.Filesystem

if !isdefined(Base, :istextmime)
    export istextmime
    istextmime(m::@compat(Union{MIME,AbstractString})) = istext(m)
end

function primarytype(t::ANY)
    tn = t.name
    if isdefined(tn, :primary)
        return tn.primary
    else
        return tn.wrapper
    end
end

if !isdefined(Base, :Threads)
    @eval module Threads
        macro threads(expr)
            return esc(expr)
        end
        threadid() = 1
        nthreads() = 1
        export @threads, threadid, nthreads
    end
    export Threads
end

if !isdefined(Base, :normalize)
    function normalize!(v::AbstractVector, p::Real=2)
        nrm = norm(v, p)
        __normalize!(v, nrm)
    end

    @inline function __normalize!(v::AbstractVector, nrm::AbstractFloat)
        #The largest positive floating point number whose inverse is less than
        #infinity
        δ = inv(prevfloat(typemax(nrm)))
        if nrm ≥ δ #Safe to multiply with inverse
            invnrm = inv(nrm)
            scale!(v, invnrm)
        else # scale elements to avoid overflow
            εδ = eps(one(nrm))/δ
            scale!(v, εδ)
            scale!(v, inv(nrm*εδ))
        end
        v
    end

    copy_oftype{T,N}(A::AbstractArray{T,N}, ::Type{T}) = copy(A)
    copy_oftype{T,N,S}(A::AbstractArray{T,N}, ::Type{S}) = convert(AbstractArray{S,N}, A)

    function normalize(v::AbstractVector, p::Real = 2)
        nrm = norm(v, p)
        if !isempty(v)
            vv = copy_oftype(v, typeof(v[1]/nrm))
            return __normalize!(vv, nrm)
        else
            T = typeof(zero(eltype(v))/nrm)
            return T[]
        end
    end

    export normalize, normalize!
end

import Base.AsyncCondition
import Base: srand, rand, rand!

if !isdefined(Base, :pointer_to_string)

    function pointer_to_string(p::Ptr{UInt8}, len::Integer, own::Bool=false)
        a = ccall(:jl_ptr_to_array_1d, Vector{UInt8},
                  (Any, Ptr{UInt8}, Csize_t, Cint), Vector{UInt8}, p, len, own)
        ccall(:jl_array_to_string, Ref{String}, (Any,), a)
    end

    pointer_to_string(p::Ptr{UInt8}, own::Bool=false) =
        pointer_to_string(p, ccall(:strlen, Csize_t, (Cstring,), p), own)

end

import Base.promote_eltype_op

import Base.LinAlg.BLAS.@blasfunc

import Base: redirect_stdin, redirect_stdout, redirect_stderr
if VERSION < v"0.6.0-dev.374"
    for (F,S) in ((:redirect_stdin, :STDIN), (:redirect_stdout, :STDOUT), (:redirect_stderr, :STDERR))
        @eval function $F(f::Function, stream)
            STDOLD = $S
            $F(stream)
            try f() finally $F(STDOLD) end
        end
    end
end

if VERSION < v"0.6.0-dev.528"
    macro __DIR__()
        Base.source_dir()
    end
    export @__DIR__
end

# PR #17302
# Provide a non-deprecated version of `@vectorize_(1|2)arg` macro which defines
# deprecated version of the function so that the depwarns can be fixed without
# breaking users.
# Packages are expected to use this to maintain the old API until all users
# of the deprecated vectorized function have migrated.
# These macros should raise a depwarn when the `0.5` support is dropped from
# `Compat` and be dropped when the support for `0.6` is dropped from `Compat`.
# Modified based on the version copied from 0.6 Base.
macro dep_vectorize_1arg(S, f)
    S = esc(S)
    f = esc(f)
    T = esc(:T)
    x = esc(:x)
    AbsArr = esc(:AbstractArray)
    ## Depwarn to be enabled when 0.5 support is dropped.
    # depwarn("Implicit vectorized function is deprecated in favor of compact broadcast syntax.",
    #         Symbol("@dep_vectorize_1arg"))
    :(@deprecate $f{$T<:$S}($x::$AbsArr{$T}) @compat($f.($x)))
end

macro dep_vectorize_2arg(S, f)
    S = esc(S)
    f = esc(f)
    T1 = esc(:T1)
    T2 = esc(:T2)
    x = esc(:x)
    y = esc(:y)
    AbsArr = esc(:AbstractArray)
    ## Depwarn to be enabled when 0.5 support is dropped.
    # depwarn("Implicit vectorized function is deprecated in favor of compact broadcast syntax.",
    #         Symbol("@dep_vectorize_2arg"))
    quote
        @deprecate $f{$T1<:$S}($x::$S, $y::$AbsArr{$T1}) @compat($f.($x,$y))
        @deprecate $f{$T1<:$S}($x::$AbsArr{$T1}, $y::$S) @compat($f.($x,$y))
        @deprecate $f{$T1<:$S,$T2<:$S}($x::$AbsArr{$T1}, $y::$AbsArr{$T2}) @compat($f.($x,$y))
    end
end

# broadcast over same length tuples, from julia#16986
if VERSION < v"0.6.0-dev.693"
    Base.Broadcast.broadcast{N}(f, t::NTuple{N}, ts::Vararg{NTuple{N}}) = map(f, t, ts...)
end

# julia#18510
if VERSION < v"0.6.0-dev.826"
    _Nullable_field2(x) = !x
else
    _Nullable_field2(x) = x
end

# julia#18484
if VERSION < v"0.6.0-dev.848"
    unsafe_get(x::Nullable) = x.value
    unsafe_get(x) = x
    export unsafe_get
    Base.isnull(x) = false
end

# julia#18977
if !isdefined(Base, :xor)
    const xor = $
    const ⊻ = xor
    export xor, ⊻
end

# julia#19246
if !isdefined(Base, :numerator)
    const numerator = num
    const denominator = den
    export numerator, denominator
end

# julia #19950
if !isdefined(Base, :iszero)
    iszero(x) = x == zero(x)
    iszero(x::Number) = x == 0
    iszero(x::AbstractArray) = all(iszero, x)
    export iszero
end

# julia　#20407
if !isdefined(Base, :(>:))
    const >: = let
        _issupertype(a::ANY, b::ANY) = issubtype(b, a)
    end
    export >:
end

# julia#19088
if VERSION < v"0.6.0-dev.1256"
    Base.take!(io::Base.AbstractIOBuffer) = takebuf_array(io)
end

# julia #17155 function composition and negation
if VERSION < v"0.6.0-dev.1883"
    export ∘
    ∘(f, g) = (x...)->f(g(x...))
    @compat Base.:!(f::Function) = (x...)->!f(x...)
end

if VERSION < v"0.6.0-dev.1632"
    # To work around unsupported syntax on Julia 0.4
    include_string("export .&, .|")
    include_string(".&(xs...) = broadcast(&, xs...)")
    include_string(".|(xs...) = broadcast(|, xs...)")
end

if VERSION < v"0.6.0-dev.2093" # Compat.isapprox to allow for NaNs
    using Base.rtoldefault
    function isapprox(x::Number, y::Number; rtol::Real=rtoldefault(x,y), atol::Real=0, nans::Bool=false)
        x == y || (isfinite(x) && isfinite(y) && abs(x-y) <= atol + rtol*max(abs(x), abs(y))) || (nans && isnan(x) && isnan(y))
    end
else
    import Base.isapprox
end

module TypeUtils
    @static if isdefined(Base, :isabstract)
        using Base: isabstract, parameter_upper_bound, typename
    else
        isabstract(t::DataType) = t.abstract
        if isdefined(Base, :TypeConstructor)
            isabstract(t::TypeConstructor) = isabstract(t.body)
        end
        isabstract(t::ANY) = false
        parameter_upper_bound(t::DataType, idx) = t.parameters[idx].ub
        typename(t::DataType) = t.name
    end
    export isabstract, parameter_upper_bound, typename
end # module TypeUtils

# @view, @views, @__dot__
include("arraymacros.jl")

# julia #18839
if VERSION < v"0.6.0-dev.1024"
    @eval module Iterators
        export countfrom, cycle, drop, enumerate, flatten, product, repeated,
               rest, take, zip, partition

        import Base: eltype, start, next, done, length, size, ndims
        using Base: tuple_type_cons
        using Base: countfrom, cycle, drop, enumerate, repeated, rest, take,
                    zip
        using Compat

        using Base: flatten
        using Base: product
        using Base: partition
    end
else
    using Base: Iterators
end

if VERSION < v"0.6.0-dev.2840"
    export IndexStyle, IndexLinear, IndexCartesian
    eval(Expr(:typealias, :IndexStyle, :(Base.LinearIndexing)))
    eval(Expr(:typealias, :IndexLinear, :(Base.LinearFast)))
    eval(Expr(:typealias, :IndexCartesian, :(Base.LinearSlow)))
    IndexStyle{T}(::Type{T}) = Base.linearindexing(T)
    IndexStyle(args...) = Base.linearindexing(args...)
end

if VERSION < v"0.6.0-dev.1653"
    for (fname, felt) in ((:zeros,:zero), (:ones,:one))
        @eval begin
            # allow signature of similar
            Base.$fname(a::AbstractArray, T::Type, dims::Tuple) = fill!(similar(a, T, dims), $felt(T))
            Base.$fname(a::AbstractArray, T::Type, dims...) = fill!(similar(a,T,dims...), $felt(T))
            Base.$fname(a::AbstractArray, T::Type=eltype(a)) = fill!(similar(a,T), $felt(T))
        end
    end
end

# https://github.com/JuliaLang/julia/pull/20203
if VERSION < v"0.6.0-dev.2283"
    # not exported
    function readline(s::IO=STDIN; chomp::Bool=true)
        if chomp
            Base.chomp!(Base.readline(s))
        else
            Base.readline(s)
        end
    end
end

# https://github.com/JuliaLang/julia/pull/18727
if VERSION < v"0.6.0-dev.838"
    Base.convert{T}(::Type{Set{T}}, s::Set{T}) = s
    Base.convert{T}(::Type{Set{T}}, s::Set) = Set{T}(s)
end

# https://github.com/JuliaLang/julia/pull/18082
if VERSION < v"0.6.0-dev.2347"
    Base.isassigned(x::Base.RefValue) = isdefined(x, :x)
end

if VERSION < v"0.6.0-dev.735"
    Base.unsafe_trunc{T<:Integer}(::Type{T}, x::Integer) = rem(x, T)
end

# https://github.com/JuliaLang/julia/pull/21346
if VERSION < v"0.6.0-pre.beta.102"
    Base.bswap(z::Complex) = Complex(bswap(real(z)), bswap(imag(z)))
end

# https://github.com/JuliaLang/julia/pull/19449
if VERSION < v"0.6.0-dev.1988"
    StringVector(n::Integer) = Vector{UInt8}(n)
else
    using Base: StringVector
end

# https://github.com/JuliaLang/julia/pull/22064
if !isdefined(Base, Symbol("@__MODULE__"))
    export @__MODULE__
    macro __MODULE__()
        return current_module()
    end
    Base.expand(mod::Module, x::ANY) = eval(mod, :(expand($(QuoteNode(x)))))
    Base.macroexpand(mod::Module, x::ANY) = eval(mod, :(macroexpand($(QuoteNode(x)))))
    Base.include_string(mod::Module, code::String, fname::String) =
        eval(mod, :(include_string($code, $fname)))
    Base.include_string(mod::Module, code::AbstractString, fname::AbstractString="string") =
        eval(mod, :(include_string($code, $fname)))
end

# https://github.com/JuliaLang/julia/pull/19784
if isdefined(Base, :invokelatest)
    import Base.invokelatest
else
    invokelatest(f, args...) = eval(current_module(), Expr(:call, f, map(QuoteNode, args)...))
end

# https://github.com/JuliaLang/julia/pull/21257
if VERSION < v"0.6.0-pre.beta.28"
    collect(A) = collect_indices(indices(A), A)
    collect_indices(::Tuple{}, A) = copy!(Array{eltype(A)}(), A)
    collect_indices(indsA::Tuple{Vararg{Base.OneTo}}, A) =
        copy!(Array{eltype(A)}(map(length, indsA)), A)
    function collect_indices(indsA, A)
        B = Array{eltype(A)}(map(length, indsA))
        copy!(B, CartesianRange(indices(B)), A, CartesianRange(indsA))
    end
else
    const collect = Base.collect
end

# https://github.com/JuliaLang/julia/pull/21378
if VERSION < v"0.6.0-pre.beta.455"
    import Base: ==, isless

    ==(x::Dates.Period, y::Dates.Period) = (==)(promote(x, y)...)
    isless(x::Dates.Period, y::Dates.Period) = isless(promote(x,y)...)

    # disallow comparing fixed to other periods
    ==(x::Dates.FixedPeriod, y::Dates.OtherPeriod) = throw(MethodError(==, (x, y)))
    ==(x::Dates.OtherPeriod, y::Dates.FixedPeriod) = throw(MethodError(==, (x, y)))
    isless(x::Dates.FixedPeriod, y::Dates.OtherPeriod) = throw(MethodError(isless, (x, y)))
    isless(x::Dates.OtherPeriod, y::Dates.FixedPeriod) = throw(MethodError(isless, (x, y)))
end

include("deprecated.jl")

# https://github.com/JuliaLang/julia/pull/21746
const macros_have_sourceloc = VERSION >= v"0.7-" && length(:(@test).args) == 2

end # module Compat


text:
using Flux
using Flux.Data: DataLoader
using Flux.Optimise: Optimiser, WeightDecay
using Flux: onehotbatch, onecold, glorot_normal, label_smoothing
using Flux.Losses: logitcrossentropy
using Statistics, Random
using Logging: with_logger
using ProgressMeter: @showprogress
import MLDatasets
using CUDA
using Formatting

using DeepUncertainty
include("utils.jl")

# function LeNet5(args; imgsize=(28, 28, 1), nclasses=10)
#     out_conv_size = (imgsize[1] ÷ 4 - 3, imgsize[2] ÷ 4 - 3, 16)

#     return Chain(
#         MCConv((5, 5), imgsize[end] => 6, args.dropout, relu),
#         MaxPool((2, 2)),
#         MCConv((5, 5), 6 => 16, args.dropout, relu),
#         MaxPool((2, 2)),
#         flatten,
#         MCDense(prod(out_conv_size), 120, args.dropout, relu),
#         MCDense(120, 84, args.dropout, relu),
#         MCDense(84, nclasses, args.dropout),
#     )
# end

function LeNet5(args; imgsize = (28, 28, 1), nclasses = 10)
    out_conv_size = (imgsize[1] ÷ 4 - 3, imgsize[2] ÷ 4 - 3, 16)
    dropout = (x; k...) -> Flux.dropout(x, args.dropout; k...)

    return Chain(
        MCLayer(Conv((5, 5), imgsize[end] => 6, relu), dropout),
        MaxPool((2, 2)),
        MCLayer(Conv((5, 5), 6 => 16, relu), dropout),
        MaxPool((2, 2)),
        flatten,
        MCLayer(Dense(prod(out_conv_size), 120, relu), dropout),
        MCLayer(Dense(120, 84, relu), dropout),
        MCLayer(Dense(84, nclasses), dropout),
    )
end

# arguments for the `train` function
Base.@kwdef mutable struct Args
    η = 3e-4             # learning rate
    λ = 0                # L2 regularizer param, implemented as weight decay
    batchsize = 256      # batch size
    epochs = 100          # number of epochs
    use_cuda = true      # if true use cuda (if available)
    infotime = 1        # report every `infotime` epochs
    checktime = 5        # Save the model every `checktime` epochs. Set to 0 for no checkpoints.
    dropout = 0.2
    sample_size = 25
end

function train(; kws...)
    args = Args(; kws...)
    use_cuda = args.use_cuda && CUDA.functional()

    if use_cuda
        device = gpu
        @info "Training on GPU"
    else
        device = cpu
        @info "Training on CPU"
    end

    ## DATA
    train_loader, test_loader, ood_train_loader, ood_test_loader = get_data(args)
    @info "Dataset MNIST: $(train_loader.nobs) train and $(test_loader.nobs) test examples"

    ## MODEL AND OPTIMIZER
    model = LeNet5(args) |> device
    @info "LeNet5 model: $(num_params(model)) trainable params"

    ps = Flux.params(model)

    opt = Nesterov(args.η)
    if args.λ > 0 # add weight decay, equivalent to L2 regularization
        opt = Optimiser(WeightDecay(args.λ), opt)
    end

    function report(epoch)
        @info "Test metrics"
        monte_carlo_evaluation(args, test_loader, model, device)
        @info "Test OOD metrics"
        monte_carlo_evaluation(args, ood_test_loader, model, device)
    end

    ## TRAINING
    @info "Start Training"
    report(0)
    for epoch = 1:args.epochs
        @showprogress for (x, y) in train_loader
            x, y = x |> device, y |> device
            gs = Flux.gradient(ps) do
                ŷ = model(x)
                loss(ŷ, y)
            end

            Flux.Optimise.update!(opt, ps, gs)
        end

        ## Printing and logging
        epoch % args.infotime == 0 && report(epoch)
    end
end

train()

ext: jl


text:
using JSON

if (length(ARGS) != 3)
    error("Arguments must have `filepath`, `year` and `sem`")
end
filepath = ARGS[1]
year = ARGS[2]
sem = ARGS[3]

stripp(x) = strip.(x)
custom_split_delim(x) = split(x, ',')
lines = readlines(filepath)
formatted_lines = stripp.(custom_split_delim.(replace.(lines, "\"" => "")))

mutable struct Course
    code::String
    name::String
    EX::Int
    A::Int
    B::Int
    C::Int
    D::Int
    P::Int
    F::Int
    X::Int
    Y::Int
    N::Int
end

Course(code, name) = Course(code, name, zeros(Int, 10, 1)...)

grades = ["A", "B", "C", "D", "P", "F", "X", "N", "Y", "EX"]

function show_course(c::Course)
    println("CourseCode: $(c.code) ($(c.name)) has grades")
    for g in grades
        println(" $(g) has grade: ", getproperty(c, Meta.parse(eval(g))))
    end
end

function format_grades_yearly_format(c::Course)
    _dict = Dict{String, Any}()
    for g in grades
        _dict[g] = getproperty(c, Meta.parse(eval(g)))
    end
    return _dict
end

courses_found = Array{Course, 1}()
cache = Dict{String, Int}()

for j in eachindex(formatted_lines)
    line = formatted_lines[j]
    print(j)
    for i in eachindex(line)
        entry = strip(line[i])
        if match(r"^[A-Z]{2}[0-9]{5}$", entry) != nothing
            if get(cache, entry, -1) == -1
                cache[entry] = length(courses_found) + 1
                push!(courses_found, Course(entry, strip(line[i+1])))
            end
            if line[i+2] ∉ grades && line[i+2] != ""
#            if entry ∈ ("ME31007",  "CY41005", "IM21001", "ME39007")
                grade = Meta.parse(eval(strip(line[i+3])))
            else
                grade = Meta.parse(eval(strip(line[i+2])))
            end
            grade == nothing && continue
            setproperty!(courses_found[cache[entry]], grade, getproperty(courses_found[cache[entry]], grade) + 1)
        end
    end
    println("")
end

all_data = JSON.parsefile("yearWiseGrades.json")
year_sem = year * sem

for c in courses_found
    if get(all_data, c.code, -1) != -1
        all_data[c.code][year_sem] = format_grades_yearly_format(c)
    else
        all_data[c.code] = Dict{String,Any}(year_sem => format_grades_yearly_format(c))
    end
end

open("yearWiseGrades.json", "w") do f
    JSON.print(f, all_data)
end

courses_info = JSON.parsefile("../courses.json")

for c in courses_found
    if get(courses_info, c.code, -1) == -1
        courses_info[c.code] = Dict{String, Any}("id" => c.code, "name" => c.name)
    end
end

open("../courses.json", "w") do f
    JSON.print(f, courses_info)
end

max_stars_repo_name: cipherLord/Kronos

max_issues_repo_name: cipherLord/Kronos

max_forks_repo_name: cipherLord/Kronos


text:
module Rules
    using ..TermInterface
    using ..Util
    import ..closure_generator

    include("patterns.jl")
    include("patterns_syntax.jl")
    include("rule_types.jl")
    include("rule_cache.jl")
    include("rule_dsl.jl")
    include("exports.jl")
end


text:
#
# Velocity gradients
#

# Diagonal
∂x_u(i, j, k, grid, u) = ∂x_caa(i, j, k, grid, u)
∂y_v(i, j, k, grid, v) = ∂y_aca(i, j, k, grid, v)
∂z_w(i, j, k, grid, w) = ∂z_aac(i, j, k, grid, w)

# Off-diagonal
∂x_v(i, j, k, grid, v) = ∂x_faa(i, j, k, grid, v)
∂x_w(i, j, k, grid, w) = ∂x_faa(i, j, k, grid, w)

∂y_u(i, j, k, grid, u) = ∂y_afa(i, j, k, grid, u)
∂y_w(i, j, k, grid, w) = ∂y_afa(i, j, k, grid, w)

∂z_u(i, j, k, grid, u) = ∂z_aaf(i, j, k, grid, u)
∂z_v(i, j, k, grid, v) = ∂z_aaf(i, j, k, grid, v)

#
# Strain components
#

# ccc strain components
Σ₁₁(i, j, k, grid, u) = ∂x_u(i, j, k, grid, u)
Σ₂₂(i, j, k, grid, v) = ∂y_v(i, j, k, grid, v)
Σ₃₃(i, j, k, grid, w) = ∂z_w(i, j, k, grid, w)

tr_Σ(i, j, k, grid, u, v, w) =
    Σ₁₁(i, j, k, grid, u) + Σ₂₂(i, j, k, grid, v) + Σ₃₃(i, j, k, grid, w)

# ffc
Σ₁₂(i, j, k, grid::Grid{T}, u, v) where T =
    T(0.5) * (∂y_u(i, j, k, grid, u) + ∂x_v(i, j, k, grid, v))

# fcf
Σ₁₃(i, j, k, grid::Grid{T}, u, w) where T =
    T(0.5) * (∂z_u(i, j, k, grid, u) + ∂x_w(i, j, k, grid, w))

# cff
Σ₂₃(i, j, k, grid::Grid{T}, v, w) where T =
    T(0.5) * (∂z_v(i, j, k, grid, v) + ∂y_w(i, j, k, grid, w))

Σ₁₂²(i, j, k, grid, u, v) = Σ₁₂(i, j, k, grid, u, v)^2
Σ₁₃²(i, j, k, grid, u, w) = Σ₁₃(i, j, k, grid, u, w)^2
Σ₂₃²(i, j, k, grid, v, w) = Σ₂₃(i, j, k, grid, v, w)^2

#
# Renamed functions for consistent function signatures
#

∂x_u(i, j, k, grid, u, v, w) = ∂x_u(i, j, k, grid, u)
∂x_v(i, j, k, grid, u, v, w) = ∂x_v(i, j, k, grid, v)
∂x_w(i, j, k, grid, u, v, w) = ∂x_w(i, j, k, grid, w)

∂y_u(i, j, k, grid, u, v, w) = ∂y_u(i, j, k, grid, u)
∂y_v(i, j, k, grid, u, v, w) = ∂y_v(i, j, k, grid, v)
∂y_w(i, j, k, grid, u, v, w) = ∂y_w(i, j, k, grid, w)

∂z_u(i, j, k, grid, u, v, w) = ∂z_u(i, j, k, grid, u)
∂z_v(i, j, k, grid, u, v, w) = ∂z_v(i, j, k, grid, v)
∂z_w(i, j, k, grid, u, v, w) = ∂z_w(i, j, k, grid, w)

Σ₁₁(i, j, k, grid, u, v, w) = Σ₁₁(i, j, k, grid, u)
Σ₂₂(i, j, k, grid, u, v, w) = Σ₂₂(i, j, k, grid, v)
Σ₃₃(i, j, k, grid, u, v, w) = Σ₃₃(i, j, k, grid, w)

Σ₁₂(i, j, k, grid, u, v, w) = Σ₁₂(i, j, k, grid, u, v)
Σ₁₃(i, j, k, grid, u, v, w) = Σ₁₃(i, j, k, grid, u, w)
Σ₂₃(i, j, k, grid, u, v, w) = Σ₂₃(i, j, k, grid, v, w)

# Symmetry relations
Σ₂₁ = Σ₁₂
Σ₃₁ = Σ₁₃
Σ₃₂ = Σ₂₃

# Trace and squared strains
tr_Σ²(ijk...) = Σ₁₁(ijk...)^2 +  Σ₂₂(ijk...)^2 +  Σ₃₃(ijk...)^2

Σ₁₂²(i, j, k, grid, u, v, w) = Σ₁₂²(i, j, k, grid, u, v)
Σ₁₃²(i, j, k, grid, u, v, w) = Σ₁₃²(i, j, k, grid, u, w)
Σ₂₃²(i, j, k, grid, u, v, w) = Σ₂₃²(i, j, k, grid, v, w)


text:
# This file is a part of JuliaFEM.
# License is MIT: see https://github.com/JuliaFEM/FEMBase.jl/blob/master/LICENSE

using FEMBase
using FEMBase.FEMTest

@testset "Poisson + Dirichlet" begin
    el1 = Element(Quad4, [1, 2, 3, 4])
    el2 = Element(Tri3, [3, 2, 5])
    el3 = Element(Seg2, [3, 5])
    el4 = Element(Seg2, [1, 4])
    X = Dict(1 => [0.0, 0.0],
             2 => [1.0, 0.0],
             3 => [1.0, 1.0],
             4 => [0.0, 1.0],
             5 => [2.0, 1.0])
    update!([el1, el2, el3, el4], "geometry", X)
    update!([el1, el2], "coefficient", 6.0)
    update!(el1, "source", 132.0)
    update!(el3, "flux", 264.0)
    update!(el4, "u 1", 0.0)
    problem = Problem(Poisson, "test problem", 1)
    bc = Problem(Dirichlet, "fixed", 1, "u")
    add_elements!(problem, [el1, el2, el3])
    add_elements!(bc, [el4])

    step = Analysis(Static)
    add_problems!(step, [problem, bc])
    ls, normu, normla = run!(step)

    @test isapprox(normu, 136.5979502042399)
    @test isapprox(normla, 280.52807346146307)

end

@testset "Poisson (including bc)" begin
    el1 = Element(Quad4, [1, 2, 3, 4])
    el2 = Element(Tri3, [3, 2, 5])
    el3 = Element(Seg2, [3, 5])
    el4 = Element(Seg2, [1, 4])
    X = Dict(1 => [0.0, 0.0],
             2 => [1.0, 0.0],
             3 => [1.0, 1.0],
             4 => [0.0, 1.0],
             5 => [2.0, 1.0])
    elements = [el1, el2, el3, el4]
    update!(elements, "geometry", X)
    update!([el1, el2], "coefficient", 6.0)
    update!(el1, "source", 132.0)
    update!(el3, "flux", 264.0)
    update!(el4, "fixed u", 0.0)
    problem = Problem(Poisson, "test problem", 1)
    add_elements!(problem, elements)

    step = Analysis(Static)
    add_problems!(step, [problem])
    ls, normu, normla = run!(step)

    @test isapprox(normu, 136.5979502042399)
    @test isapprox(normla, 569.2099788303083)
end

@testset "Test @test_resource" begin
    fn = @test_resource("mesh.inp")
    @test occursin("mesh.inp", fn)
end

@testset "Test read_mtx_from_file and read_mtx_from_string" begin
    data = """
    1,1, 1,1,  1.0
    """
    K = read_mtx_from_string(data)
    @test isapprox(K[1,1], 1.0)
    K2 = read_mtx_from_file(@test_resource("eb3arxs4_STIF1.mtx"))
    @test isapprox(K2[1,1], 1.8e6)
end


text:
using Pkg
Pkg.activate(".")

using Genie, HTTP
import Genie.Router: route, POST, @params
import Genie.Requests: jsonpayload

route("/json-error", method = POST) do
  error("500, sorry")
end

Genie.AppServer.startup()

@test_throws HTTP.ExceptionRequest.StatusError HTTP.request("POST", "http://localhost:8000/json-error", [("Content-Type", "application/json; charset=utf-8")], """{"greeting":"hello"}""")

exit(0)


text:
# This file is a part of Pukeko.jl.
# License is MIT. https://github.com/IainNZ/Pukeko.jl
# This file is adapted from the Julia `test/int.jl` file.
# It serves as an example.
# License is MIT: https://julialang.org/license

module BaseInt
    using Pukeko

    function test_flipsign_copysign()
        for y in (-4, Float32(-4), -4.0, big(-4.0))
            @test flipsign(3, y)  == -3
            @test flipsign(-3, y) == 3
            @test copysign(3, y)  == -3
            @test copysign(-3, y) == -3
        end

        for y in (4, Float32(4), 4.0, big(4.0))
            @test flipsign(3, y)  == 3
            @test flipsign(-3, y) == -3
            @test copysign(3, y)  == 3
            @test copysign(-3, y) == 3
        end
    end

    function _flipsign_copysign_type(U)
        for T in (Base.BitInteger_types..., BigInt,
                  Rational{Int}, Rational{BigInt},
                  Float16, Float32, Float64)
            @test typeof(copysign(T(3), U(4))) === T
            @test typeof(flipsign(T(3), U(4))) === T
        end
        # Bool promotes to Int
        U <: Unsigned && return
        for x in [true, false]
            @test flipsign(x, U(4)) === Int(x)
            @test flipsign(x, U(-1)) === -Int(x)
            @test copysign(x, U(4)) === Int(x)
            @test copysign(x, U(-1)) === -Int(x)
        end
    end
    @parametric _flipsign_copysign_type (Base.BitInteger_types..., BigInt,
                                         Rational{Int}, Rational{BigInt},
                                         Float16, Float32, Float64)

    function _flipsign_copysign_typemin(T)
        for U in (Base.BitSigned_types..., BigInt, Float16, Float32, Float64)
            @test flipsign(typemin(T), U(-1)) == typemin(T)
            @test copysign(typemin(T), U(-1)) == typemin(T)
        end
    end
    @parametric _flipsign_copysign_typemin Base.BitInteger_types

    function test_flipsign_with_floats()
        for s1 in (-1,+1), s2 in (-1,+1)
            @test flipsign(Int16(3s1), Float16(3s2)) === Int16(3s1*s2)
            @test flipsign(Int32(3s1), Float32(3s2)) === Int32(3s1*s2)
            @test flipsign(Int64(3s1), Float64(3s2)) === Int64(3s1*s2)
        end
    end
end

import Pukeko
Pukeko.run_tests(BaseInt)


text:
using Random: GLOBAL_RNG

struct ELBO <: VariationalObjective end

function (elbo::ELBO)(alg, q, logπ, num_samples; kwargs...)
    return elbo(GLOBAL_RNG, alg, q, logπ, num_samples; kwargs...)
end

const elbo = ELBO()


text:
using Test, Random, StructArrays, Tables, StaticArrays, OffsetArrays
using GeometryBasics
using LinearAlgebra
using GeometryBasics: attributes

@testset "GeometryBasics" begin

@testset "algorithms" begin
    cube = Rect(Vec3f0(-0.5), Vec3f0(1))
    cube_faces = decompose(TriangleFace{Int}, faces(cube))
    cube_vertices = decompose(Point{3,Float32}, cube)
    @test area(cube_vertices, cube_faces) == 6
    mesh = Mesh(cube_vertices, cube_faces)
    @test GeometryBasics.volume(mesh) ≈ 1
    @test GeometryBasics.volume(cube) ≈ 1

    rect = Rect(1, 2, 7.5, 2.0)
    @test GeometryBasics.area(rect) ≈ 15

    points_cwise = Point2f0[(0,0), (0,1), (1,1)]
    points_ccwise = Point2f0[(0,0), (1,0), (1,1)]
    @test area(points_cwise) ≈ -0.5
    @test area(points_ccwise) ≈ 0.5
    @test area(OffsetArray(points_cwise, -2)) ≈ -0.5

    points3d = Point3f0[(0,0,0), (0,0,1), (0,1,1)]
    @test area(OffsetArray(points3d, -2)) ≈ 0.5

    pm2d = [PointMeta(0.0, 0.0, a=:d), PointMeta(0.0, 1.0, a=:e), PointMeta(1.0, 0.0, a=:f)]
    @test area(pm2d) ≈ -0.5

    pm3d = [PointMeta(0.0, 0.0, 0.0, a=:d), PointMeta(0.0, 1.0, 0.0, a=:e), PointMeta(1.0, 0.0, 0.0, a=:f)]
    @test_broken area(pm3d) ≈ 0.5 # Currently broken as zero(PointMeta(0.0, 0.0, 0.0, a=:d)) fails
end

@testset "embedding metadata" begin
    @testset "Meshes" begin

        @testset "per vertex attributes" begin
            points = rand(Point{3, Float64}, 8)
            tfaces = TetrahedronFace{Int}[(1, 2, 3, 4), (5, 6, 7, 8)]
            normals = rand(SVector{3, Float64}, 8)
            stress = LinRange(0, 1, 8)
            mesh = Mesh(meta(points, normals = normals, stress = stress), tfaces)

            @test hasproperty(coordinates(mesh), :stress)
            @test hasproperty(coordinates(mesh), :normals)
            @test coordinates(mesh).stress === stress
            @test coordinates(mesh).normals === normals
            @test coordinates(mesh).normals === normals
            @test GeometryBasics.faces(mesh) === tfaces
            @test propertynames(coordinates(mesh)) == (:position, :normals, :stress)

        end

        @testset "per face attributes" begin

            # Construct a cube out of Quads
            points = Point{3, Float64}[
                (0.0, 0.0, 0.0), (2.0, 0.0, 0.0),
                (2.0, 2.0, 0.0), (0.0, 2.0, 0.0),
                (0.0, 0.0, 12.0), (2.0, 0.0, 12.0),
                (2.0, 2.0, 12.0), (0.0, 2.0, 12.0)
            ]

            facets = QuadFace{Cint}[
                1:4,
                5:8,
                [1,5,6,2],
                [2,6,7,3],
                [3, 7, 8, 4],
                [4, 8, 5, 1]
            ]

            markers = Cint[-1, -2, 0, 0, 0, 0]
            # attach some additional information to our faces!
            mesh = Mesh(points, meta(facets, markers = markers))
            @test hasproperty(GeometryBasics.faces(mesh), :markers)
            # test with === to assert we're not doing any copies
            @test GeometryBasics.faces(mesh).markers === markers
            @test coordinates(mesh) === points
            @test metafree(GeometryBasics.faces(mesh)) === facets

        end

    end

    @testset "polygon with metadata" begin
        polys = [Polygon(rand(Point{2, Float32}, 20)) for i in 1:10]
        pnames = [randstring(4) for i in 1:10]
        numbers = LinRange(0.0, 1.0, 10)
        bin = rand(Bool, 10)
        # create a polygon
        poly = PolygonMeta(polys[1], name = pnames[1], value = numbers[1], category = bin[1])
        # create a MultiPolygon with the right type & meta information!
        multipoly = MultiPolygonMeta(polys, name = pnames, value = numbers, category = bin)
        @test multipoly isa AbstractVector
        @test poly isa GeometryBasics.AbstractPolygon

        @test GeometryBasics.getcolumn(poly, :name) == pnames[1]
        @test GeometryBasics.MetaFree(PolygonMeta) == Polygon

        @test GeometryBasics.getcolumn(multipoly, :name) == pnames
        @test GeometryBasics.MetaFree(MultiPolygonMeta) == MultiPolygon

        meta_p = meta(polys[1], boundingbox=Rect(0, 0, 2, 2))
        @test meta_p.boundingbox === Rect(0, 0, 2, 2)
        @test metafree(meta_p) === polys[1]
        attributes(meta_p) == Dict{Symbol, Any}(:boundingbox => meta_p.boundingbox,
                                                :polygon => polys[1])
    end
    @testset "point with metadata" begin
        p = Point(1.1, 2.2)
        @test p isa AbstractVector{Float64}
        pm = PointMeta(1.1, 2.2; a=1, b=2)
        p1 = Point(2.2, 3.6)
        p2 = [p, p1]
        @test coordinates(p2) == p2
        @test meta(pm) === (a=1, b=2)
        @test metafree(pm) === p
        @test propertynames(pm) == (:position, :a, :b)
        @test GeometryBasics.MetaFree(typeof(pm)) == Point{2,Float64}
        @test_broken zero(pm) == [0, 0]
    end

    @testset "MultiPoint with metadata" begin
        p = collect(Point{2, Float64}(x, x+1) for x in 1:5)
        @test p isa AbstractVector
        mpm = MultiPointMeta(p, a=1, b=2)
        @test coordinates(mpm) == mpm
        @test meta(mpm) === (a=1, b=2)
        @test metafree(mpm) == p
        @test propertynames(mpm) == (:points, :a, :b)
    end

    @testset "LineString with metadata" begin
        linestring = LineStringMeta(Point{2, Int}[(10, 10), (20, 20), (10, 40)], a = 1, b = 2)
        @test linestring isa AbstractVector
        @test meta(linestring) === (a = 1, b = 2)
        @test metafree(linestring) == linestring
        @test propertynames(linestring) == (:lines, :a, :b)
    end

    @testset "MultiLineString with metadata" begin
        linestring1 = LineString(Point{2, Int}[(10, 10), (20, 20), (10, 40)])
        linestring2 = LineString(Point{2, Int}[(40, 40), (30, 30), (40, 20), (30, 10)])
        multilinestring = MultiLineString([linestring1, linestring2])
        multilinestringmeta = MultiLineStringMeta([linestring1, linestring2]; boundingbox = Rect(1.0, 1.0, 2.0, 2.0))
        @test multilinestringmeta isa AbstractVector
        @test meta(multilinestringmeta) === (boundingbox = Rect(1.0, 1.0, 2.0, 2.0),)
        @test metafree(multilinestringmeta) == multilinestring
        @test propertynames(multilinestringmeta) == (:linestrings, :boundingbox)
    end

    @testset "Mesh with metadata" begin
       m = triangle_mesh(Sphere(Point3f0(0), 1))
       m_meta = MeshMeta(m; boundingbox=Rect(1.0, 1.0, 2.0, 2.0))
       @test meta(m_meta) === (boundingbox = Rect(1.0, 1.0, 2.0, 2.0),)
       @test metafree(m_meta) === m
       @test propertynames(m_meta) == (:mesh, :boundingbox)
   end
end

@testset "embedding MetaT" begin
    @testset "MetaT{Polygon}" begin
        polys = [Polygon(rand(Point{2, Float32}, 20)) for i in 1:10]
        multipol = MultiPolygon(polys)
        pnames = [randstring(4) for i in 1:10]
        numbers = LinRange(0.0, 1.0, 10)
        bin = rand(Bool, 10)
        # create a polygon
        poly = MetaT(polys[1], name = pnames[1], value = numbers[1], category = bin[1])
        # create a MultiPolygon with the right type & meta information!
        multipoly = MetaT(multipol, name = pnames, value = numbers, category = bin)
        @test multipoly isa MetaT
        @test poly isa MetaT

        @test GeometryBasics.getcolumn(poly, :name) == pnames[1]
        @test GeometryBasics.getcolumn(multipoly, :name) == pnames

        meta_p = MetaT(polys[1], boundingbox=Rect(0, 0, 2, 2))
        @test meta_p.boundingbox === Rect(0, 0, 2, 2)
        @test GeometryBasics.metafree(meta_p) == polys[1]
        @test GeometryBasics.metafree(poly) == polys[1]
        @test GeometryBasics.metafree(multipoly) == multipol
        @test GeometryBasics.meta(meta_p) == (boundingbox = GeometryBasics.HyperRectangle{2,Int64}([0, 0], [2, 2]),)
        @test GeometryBasics.meta(poly) == (name = pnames[1], value = 0.0, category = bin[1])
        @test GeometryBasics.meta(multipoly) == (name = pnames, value = numbers, category = bin)
    end

    @testset "MetaT{Point}" begin
        p = Point(1.1, 2.2)
        @test p isa AbstractVector{Float64}
        pm = MetaT(Point(1.1, 2.2); a=1, b=2)
        p1 = Point(2.2, 3.6)
        p2 = [p, p1]
        @test coordinates(p2) == p2
        @test pm.meta === (a=1, b=2)
        @test pm.main === p
        @test propertynames(pm) == (:main, :a, :b)
        @test GeometryBasics.metafree(pm) == p
        @test GeometryBasics.meta(pm) == (a = 1, b = 2)
    end

    @testset "MetaT{MultiPoint}" begin
        p = collect(Point{2, Float64}(x, x+1) for x in 1:5)
        @test p isa AbstractVector
        mpm = MetaT(MultiPoint(p); a=1, b=2)
        @test coordinates(mpm.main) == Point{2, Float64}[(x, x+1) for x in 1:5]
        @test mpm.meta === (a=1, b=2)
        @test mpm.main == p
        @test propertynames(mpm) == (:main, :a, :b)
        @test GeometryBasics.metafree(mpm) == p
        @test GeometryBasics.meta(mpm) == (a = 1, b = 2)
    end

    @testset "MetaT{LineString}" begin
        linestring = MetaT(LineString(Point{2, Int}[(10, 10), (20, 20), (10, 40)]), a = 1, b = 2)
        @test linestring isa MetaT
        @test linestring.meta === (a = 1, b = 2)
        @test propertynames(linestring) == (:main, :a, :b)
        @test GeometryBasics.metafree(linestring) == LineString(Point{2, Int}[(10, 10), (20, 20), (10, 40)])
        @test GeometryBasics.meta(linestring) == (a = 1, b = 2)
    end

    @testset "MetaT{MultiLineString}" begin
        linestring1 = LineString(Point{2, Int}[(10, 10), (20, 20), (10, 40)])
        linestring2 = LineString(Point{2, Int}[(40, 40), (30, 30), (40, 20), (30, 10)])
        multilinestring = MultiLineString([linestring1, linestring2])
        multilinestringmeta = MetaT(MultiLineString([linestring1, linestring2]); boundingbox = Rect(1.0, 1.0, 2.0, 2.0))
        @test multilinestringmeta isa MetaT
        @test multilinestringmeta.meta === (boundingbox = Rect(1.0, 1.0, 2.0, 2.0),)
        @test multilinestringmeta.main == multilinestring
        @test propertynames(multilinestringmeta) == (:main, :boundingbox)
        @test GeometryBasics.metafree(multilinestringmeta) == multilinestring
        @test GeometryBasics.meta(multilinestringmeta) == (boundingbox = GeometryBasics.HyperRectangle{2,Float64}([1.0, 1.0], [2.0, 2.0]),)
    end

    #=
    So mesh works differently for MetaT
    since `MetaT{Point}` not subtyped to `AbstractPoint`
    =#

   @testset "MetaT{Mesh}" begin
        @testset "per vertex attributes" begin
            points = rand(Point{3, Float64}, 8)
            tfaces = TetrahedronFace{Int}[(1, 2, 3, 4), (5, 6, 7, 8)]
            normals = rand(SVector{3, Float64}, 8)
            stress = LinRange(0, 1, 8)
            mesh_nometa = Mesh(points, tfaces)
            mesh = MetaT(mesh_nometa, normals = normals, stress = stress)

            @test hasproperty(mesh, :stress)
            @test hasproperty(mesh, :normals)
            @test mesh.stress == stress
            @test mesh.normals == normals
            @test GeometryBasics.faces(mesh.main) == tfaces
            @test propertynames(mesh) == (:main, :normals, :stress)
        end
    end
end

@testset "view" begin
    @testset "TupleView" begin
        x = [1, 2, 3, 4, 5, 6]
        y = TupleView{2, 1}(x)
        @test y == [(1, 2), (2, 3), (3, 4), (4, 5), (5, 6)]

        y = TupleView{2}(x)
        @test y == [(1, 2), (3, 4), (5, 6)]

        y = TupleView{2, 3}(x)
        @test y == [(1, 2), (4, 5)]

        y = TupleView{3, 1}(x)
        @test y == [(1, 2, 3), (2, 3, 4), (3, 4, 5), (4, 5, 6)]

        y = TupleView{2, 1}(x, connect = true)
        @test y == [(1, 2), (2, 3), (3, 4), (4, 5), (5, 6), (6, 1)]

    end

    @testset "connected views" begin
        numbers = [1, 2, 3, 4, 5, 6]
        x = connect(numbers, Point{2})

        @test x == Point[(1, 2), (3, 4), (5, 6)]

        line = connect(x, Line, 1)
        @test line == [Line(Point(1, 2), Point(3, 4)), Line(Point(3, 4), Point(5, 6))]

        triangles = connect(x, Triangle)
        @test triangles == [Triangle(Point(1, 2), Point(3, 4), Point(5, 6))]
        x = connect([1, 2, 3, 4, 5, 6, 7, 8], Point{2})
        tetrahedra = connect(x, NSimplex{4})
        @test tetrahedra == [Tetrahedron(x[1], x[2], x[3], x[4])]

        @testset "matrix non-copy point views" begin
            # point in row
            points = [1 2; 1 4; 66 77]
            comparison = [Point(1, 2), Point(1, 4), Point(66, 77)]
            @test connect(points, Point{2}) == comparison
            # point in column
            points = [1 1 66; 2 4 77]
            # huh, reinterpret array doesn't seem to like `==`
            @test all(((a,b),)-> a==b, zip(connect(points, Point{2}), comparison))
        end
    end

    @testset "face views" begin
        numbers = [1, 2, 3, 4, 5, 6]
        points = connect(numbers, Point{2})
        faces = connect([1, 2, 3], TriangleFace)
        triangles = connect(points, faces)
        @test triangles == [Triangle(Point(1, 2), Point(3, 4), Point(5, 6))]
        x = Point{3}(1.0)
        triangles = connect([x], [TriangleFace(1, 1, 1)])
        @test triangles == [Triangle(x, x, x)]
        points = connect([1, 2, 3, 4, 5, 6, 7, 8], Point{2})
        faces = connect([1, 2, 3, 4], SimplexFace{4})
        triangles = connect(points, faces)
        @test triangles == [Tetrahedron(points...)]
    end

    @testset "reinterpret" begin
        numbers = collect(reshape(1:6, 2, 3))
        points = reinterpret(Point{2, Int}, numbers)
        @test points[1] === Point(1, 2)
        @test points[2] === Point(3, 4)
        numbers[4] = 0
        @test points[2] === Point(3, 0)
    end
end

@testset "constructors" begin
    @testset "LineFace" begin

        points = connect([1, 2, 3, 4, 5, 6], Point{2})
        linestring = LineString(points)
        @test linestring == [Line(points[1], points[2]), Line(points[2], points[3])]

        points = rand(Point{2, Float64}, 4)
        linestring = LineString(points, 2)
        @test linestring == [Line(points[1], points[2]), Line(points[3], points[4])]

        linestring = LineString([points[1] => points[2], points[2] => points[3]])
        @test linestring == [Line(points[1], points[2]), Line(points[2], points[3])]

        faces = [1, 2, 3]
        linestring = LineString(points, faces)
        @test linestring == LineString([points[1] => points[2], points[2] => points[3]])
        a, b, c, d = Point(1, 2), Point(3, 4), Point(5, 6), Point(7, 8)
        points = [a, b, c, d]; faces = [1, 2, 3, 4]
        linestring = LineString(points, faces, 2)
        @test linestring == LineString([a => b, c => d])

        faces = [LineFace(1, 2)
        , LineFace(3, 4)]
        linestring = LineString(points, faces)
        @test linestring == LineString([a => b, c => d])
    end

    @testset "Polygon" begin

        points = connect([1, 2, 3, 4, 5, 6], Point{2})
        polygon = Polygon(points)
        @test polygon == Polygon(LineString(points))

        points = rand(Point{2, Float64}, 4)
        linestring = LineString(points, 2)
        @test Polygon(points, 2) == Polygon(linestring)

        faces = [1, 2, 3]
        polygon = Polygon(points, faces)
        @test polygon == Polygon(LineString(points, faces))

        a, b, c, d = Point(1, 2), Point(3, 4), Point(5, 6), Point(7, 8)
        points = [a, b, c, d]; faces = [1, 2, 3, 4]
        polygon = Polygon(points, faces, 2)
        @test polygon == Polygon(LineString(points, faces, 2))

        faces = [LineFace(1, 2), LineFace(3, 4)]
        polygon = Polygon(points, faces)
        @test polygon == Polygon(LineString(points, faces))
        @test ndims(polygon) === 2
    end

    @testset "Mesh" begin

        numbers = [1, 2, 3, 4, 5, 6]
        points = connect(numbers, Point{2})

        mesh = Mesh(points, [1,2,3])
        @test mesh == [Triangle(points...)]

        x = Point{3}(1.0)
        mesh = Mesh([x], [TriangleFace(1, 1, 1)])
        @test mesh == [Triangle(x, x, x)]

        points = connect([1, 2, 3, 4, 5, 6, 7, 8], Point{2})
        faces = connect([1, 2, 3, 4], SimplexFace{4})
        mesh = Mesh(points, faces)
        @test mesh == [Tetrahedron(points...)]

        points = rand(Point3f0, 8)
        tfaces = [GLTriangleFace(1, 2, 3), GLTriangleFace(5, 6, 7)]
        normals = rand(Vec3f0, 8)
        uv = rand(Vec2f0, 8)
        mesh = Mesh(points, tfaces)
        meshuv = Mesh(meta(points; uv=uv), tfaces)
        meshuvnormal = Mesh(meta(points; normals=normals, uv=uv), tfaces)

        @test mesh isa GLPlainMesh
        @test meshuv isa GLUVMesh3D
        @test meshuvnormal isa GLNormalUVMesh3D

        t = Tesselation(FRect2D(0, 0, 2, 2), (30, 30))
        m = GeometryBasics.mesh(t, pointtype=Point3f0, facetype=QuadFace{Int})
        m2 = GeometryBasics.mesh(m, facetype=QuadFace{GLIndex})
        @test GeometryBasics.faces(m2) isa Vector{QuadFace{GLIndex}}
        @test GeometryBasics.coordinates(m2) isa Vector{Point3f0}

    end

    @testset "Multi geometries" begin
        # coordinates from https://en.wikipedia.org/wiki/Well-known_text_representation_of_geometry#Geometric_objects
        points = Point{2, Int}[(10, 40), (40, 30), (20, 20), (30, 10)]
        multipoint = MultiPoint(points)
        @test size(multipoint) === size(points)
        @test multipoint[3] === points[3]

        linestring1 = LineString(Point{2, Int}[(10, 10), (20, 20), (10, 40)])
        linestring2 = LineString(Point{2, Int}[(40, 40), (30, 30), (40, 20), (30, 10)])
        multilinestring = MultiLineString([linestring1, linestring2])
        @test size(multilinestring) === (2,)
        @test multilinestring[1] === linestring1
        @test multilinestring[2] === linestring2

        polygon11 = Polygon(Point{2, Int}[(30, 20), (45, 40), (10, 40), (30, 20)])
        polygon12 = Polygon(Point{2, Int}[(15, 5), (40, 10), (10, 20), (5, 10), (15, 5)])
        multipolygon1 = MultiPolygon([polygon11, polygon12])
        @test size(multipolygon1) === (2,)
        @test multipolygon1[1] === polygon11
        @test multipolygon1[2] === polygon12

        polygon21 = Polygon(Point{2, Int}[(40, 40), (20, 45), (45, 30), (40, 40)])
        polygon22 = Polygon(LineString(Point{2, Int}[(20, 35), (10, 30), (10, 10), (30, 5), (45, 20), (20, 35)]),
            [LineString(Point{2, Int}[(30, 20), (20, 15), (20, 25), (30, 20)])])
        multipolygon2 = MultiPolygon([polygon21, polygon22])
        @test size(multipolygon2) === (2,)
        @test multipolygon2[1] === polygon21
        @test multipolygon2[2] === polygon22
    end

end

@testset "decompose/triangulation" begin
    primitive = Sphere(Point3f0(0), 1)
    @test ndims(primitive) === 3
    mesh = triangle_mesh(primitive)
    @test decompose(Point, mesh) isa Vector{Point3f0}
    @test decompose(Point, primitive) isa Vector{Point3f0}

    primitive = Rect2D(0, 0, 1, 1)
    mesh = triangle_mesh(primitive)

    @test decompose(Point, mesh) isa Vector{Point2f0}
    @test decompose(Point, primitive) isa Vector{Point2{Int}}

    primitive = Rect3D(0, 0, 0, 1, 1, 1)
    triangle_mesh(primitive)

    primitive = Sphere(Point3f0(0), 1)
    m_normal = normal_mesh(primitive)
    @test normals(m_normal) isa Vector{Vec3f0}
    primitive = Rect2D(0, 0, 1, 1)
    m_normal = normal_mesh(primitive)
    @test normals(m_normal) isa Vector{Vec3f0}
    primitive = Rect3D(0, 0, 0, 1, 1, 1)
    m_normal = normal_mesh(primitive)
    @test normals(m_normal) isa Vector{Vec3f0}

    points = decompose(Point2f0, Circle(Point2f0(0), 1))
    tmesh = triangle_mesh(points)
    @test normals(tmesh) == nothing

    m = GeometryBasics.mesh(Sphere(Point3f0(0), 1))
    @test normals(m) == nothing
    m_normals = pointmeta(m, Normal())
    @test normals(m_normals) isa Vector{Vec3f0}

    @test texturecoordinates(m) == nothing
    r2 = Rect2D(0.0, 0.0, 1.0, 1.0)
    @test iterate(texturecoordinates(r2)) == ((0.0, 1.0), ((0.0, 2), (1.0, 2)))
    r3 = Rect3D(0.0, 0.0, 1.0, 1.0, 2.0, 2.0)
    @test iterate(texturecoordinates(r3)) == ([0, 0, 0], 2)
    uv = decompose_uv(m)
    @test Rect(Point.(uv)) == Rect(0, 0, 1, 1)

    points = decompose(Point2f0, Circle(Point2f0(0), 1))
    m = GeometryBasics.mesh(points)
    @test coordinates(m) === points

    linestring = LineString(Point{2, Int}[(10, 10), (20, 20), (10, 40)])
    pts = Point{2, Int}[(10, 10), (20, 20), (10, 40)]
    linestring = LineString(pts)
    pts_decomp = decompose(Point{2, Int}, linestring)
    @test pts === pts_decomp

    pts_ext = Point{2, Int}[(5, 1), (3, 3), (4, 8), (1, 2), (5, 1)]
    ls_ext = LineString(pts_ext)
    pts_int1 = Point{2, Int}[(2, 2), (3, 8),(5, 6), (3, 4), (2, 2)]
    ls_int1 = LineString(pts_int1)
    pts_int2 = Point{2, Int}[(3, 2), (4, 5),(6, 1), (1, 4), (3, 2)]
    ls_int2 =  LineString(pts_int2)
    poly_ext = Polygon(ls_ext)
    poly_ext_int = Polygon(ls_ext, [ls_int1, ls_int2])
    @test decompose(Point{2, Int}, poly_ext) == pts_ext
    @test decompose(Point{2, Int}, poly_ext_int) == [pts_ext..., pts_int1..., pts_int2...]
end

@testset "mesh" begin
    primitive = Triangle(Point2f0(0), Point2f0(1), Point2f0(1,0))
    m = GeometryBasics.mesh(primitive)
    @test length(faces(m)) == 1
end

@testset "convert mesh + meta" begin
    m = uv_normal_mesh(Circle(Point2f0(0), 1f0))
    # for 2D primitives we dont actually calculate normals
    @test !hasproperty(m, :normals)
end

@testset "convert mesh + meta" begin
    m = uv_normal_mesh(FRect3D(Vec3f0(-1), Vec3f0(1, 2, 3)))
    m_normal = normal_mesh(m)
    # make sure we don't loose the uv
    @test hasproperty(m_normal, :uv)
    @test m == m_normal
    # Make sure we don't create any copies
    @test m.position === m_normal.position
    @test m.normals === m_normal.normals
    @test m.uv === m_normal.uv

    m = GeometryBasics.mesh(FRect3D(Vec3f0(-1), Vec3f0(1, 2, 3));
                            uv=Vec2{Float64}, normaltype=Vec3{Float64}, pointtype=Point3{Float64})
    m_normal = normal_mesh(m)
    @test hasproperty(m_normal, :uv)
    @test m.position !== m_normal.position
    @test m.normals !== m_normal.normals
    # uv stays untouched, since we don't specify the element type in normalmesh
    @test m.uv === m_normal.uv
end

@testset "modifying meta" begin
    xx = rand(10)
    points = rand(Point3f0, 10)
    m = GeometryBasics.Mesh(meta(points, xx=xx), GLTriangleFace[(1,2,3), (3,4,5)])
    color = rand(10)
    m = pointmeta(m; color=color)

    @test hasproperty(m, :xx)
    @test hasproperty(m, :color)
    @test_throws ErrorException GeometryBasics.MetaType(Simplex)
    @test_throws ErrorException GeometryBasics.MetaFree(Simplex)

    @test m.xx === xx
    @test m.color === color

    m, colpopt = GeometryBasics.pop_pointmeta(m, :color)
    m, xxpopt = GeometryBasics.pop_pointmeta(m, :xx)

    @test propertynames(m) == (:position,)
    @test colpopt === color
    @test xxpopt === xx

    @testset "creating meta" begin
        x = Point3f0[(1,3,4)]
        # no meta gets added, so should stay the same
        @test meta(x) === x
        @test meta(x, value=[1]).position === x
    end
    pos = Point2f0[(10, 2)]
    m = Mesh(meta(pos, uv=[Vec2f0(1, 1)]), [GLTriangleFace(1, 1, 1)])
    @test m.position === pos
end

@testset "mesh conversion" begin
    s = Sphere(Point3(0.0), 1.0)
    m = GeometryBasics.mesh(s)
    @test m isa Mesh{3, Float64}
    @test coordinates(m) isa Vector{Point{3, Float64}}
    @test GeometryBasics.faces(m) isa Vector{GLTriangleFace}
    # Check, that decompose isn't making a copy for matching eltype
    @test coordinates(m) === decompose(Point{3, Float64}, m)

    tmesh = triangle_mesh(m)
    @test tmesh isa GLPlainMesh
    @test coordinates(tmesh) === decompose(Point3f0, tmesh)

    nmesh = normal_mesh(m)
    @test nmesh isa GLNormalMesh
    @test metafree(coordinates(nmesh)) === decompose(Point3f0, nmesh)
    @test normals(nmesh) === decompose_normals(nmesh)

    m = GeometryBasics.mesh(s, pointtype=Point3f0)
    @test m isa Mesh{3, Float32}
    @test coordinates(m) isa Vector{Point3f0}
    @test GeometryBasics.faces(m) isa Vector{GLTriangleFace}
end

@testset "lines intersects" begin
    a = Line(Point(0.0, 0.0), Point(4.0, 1.0))
    b = Line(Point(0.0, 0.25), Point(3.0, 0.25))
    c = Line(Point(0.0, 0.25), Point(0.5, 0.25))
    d = Line(Point(0.0, 0.0), Point(0.0, 4.0))
    e = Line(Point(1.0, 0.0), Point(0.0, 4.0))
    f = Line(Point(5.0, 0.0), Point(6.0, 0.0))

    @test intersects(a, b) === (true, Point(1.0, 0.25))
    @test intersects(a, c) === (false, Point(0.0, 0.0))
    @test intersects(d, d) === (false, Point(0.0, 0.0))
    found, point = intersects(d, e)
    @test found && point ≈ Point(0.0, 4.0)
    @test intersects(a, f) === (false, Point(0.0, 0.0))
end

@testset "Offsetintegers" begin
    x = 1
    @test GeometryBasics.raw(x) isa Int64
    @test GeometryBasics.value(x) == x

    x = ZeroIndex(1)
    @test eltype(x) == Int64

    x = OffsetInteger{0}(1)
    @test typeof(x) == OffsetInteger{0,Int64}

    x1 = OffsetInteger{0}(2)
    @test GeometryBasics.pure_max(x, x1) == x1
    @test promote_rule(typeof(x), typeof(x1)) == OffsetInteger{0,Int64}
    x2 = 1
    @test promote_rule(typeof(x2), typeof(x1)) == Int64
    @test Base.to_index(x1) == 2
    @test -(x1) == OffsetInteger{0,Int64}(-2)
    @test abs(x1) == OffsetInteger{0,Int64}(2)
    @test +(x, x1) == OffsetInteger{0,Int64}(3)
    @test *(x, x1) == OffsetInteger{0,Int64}(2)
    @test -(x, x1) == OffsetInteger{0,Int64}(-1)
    #test for /
    @test div(x, x1) == OffsetInteger{0,Int64}(0)
    @test !==(x, x1)
    @test !>=(x, x1)
    @test <=(x, x1)
    @test !>(x, x1)
    @test <(x, x1)
end

@testset "MetaT and heterogeneous data" begin
    ls = [LineString([Point(i, (i+1)^2/6), Point(i*0.86,i+5), Point(i/3, i/7)]) for i in 1:10]
    mls = MultiLineString([LineString([Point(i+1, (i)^2/6), Point(i*0.75,i+8), Point(i/2.5, i/6.79)]) for i in 5:10])
    poly = Polygon(Point{2, Int}[(40, 40), (20, 45), (45, 30), (40, 40)])
    geom = [ls..., mls, poly]
    prop = Any[(country_states = "India$(i)", rainfall = (i*9)/2) for i in 1:11]
    push!(prop, (country_states = 12, rainfall = 1000))   # a pinch of heterogeneity

    feat = [MetaT(i, j) for (i,j) = zip(geom, prop)]
    sa = meta_table(feat)

    @test nameof(eltype(feat)) == :MetaT
    @test eltype(sa) === MetaT{Any,(:country_states, :rainfall),Tuple{Any,Float64}}
    @test propertynames(sa) === (:main, :country_states, :rainfall)
    @test getproperty(sa, :country_states) isa Array{Any}
    @test getproperty(sa, :main) == geom

    maintype, metanames, metatype = GeometryBasics.getnamestypes(typeof(feat[1]))
    @test (metanames, metatype) == ((:country_states, :rainfall), Tuple{String,Float64})

    @test StructArrays.createinstance(typeof(feat[1]), LineString([Point(1, (2)^2/6), Point(1*0.86,6), Point(1/3, 1/7)]), "Mumbai", 100) isa typeof(feat[1])

    @test Base.getindex(feat[1], 1) isa Line
    @test Base.size(feat[1]) == (2,)
end

@testset "StructArrays integration" begin
    pt = meta(Point(0.0, 0.0), color="red", alpha=0.1)
    @test StructArrays.component(pt, :position) == Point(0.0, 0.0)
    @test StructArrays.component(pt, :color) == "red"
    @test StructArrays.component(pt, :alpha) == 0.1
    @test StructArrays.staticschema(typeof(pt)) ==
        NamedTuple{(:position, :color, :alpha), Tuple{Point2{Float64}, String, Float64}}
    @test StructArrays.createinstance(typeof(pt), Point(0.0, 0.0), "red", 0.1) == pt

    s = StructArray([pt, pt])
    @test StructArrays.components(s) == (
        position = [Point(0.0, 0.0), Point(0.0, 0.0)],
        color = ["red", "red"],
        alpha = [0.1, 0.1]
    )
    s[2] = meta(Point(0.1, 0.1), color="blue", alpha=0.3)
    @test StructArrays.components(s) == (
        position = [Point(0.0, 0.0), Point(0.1, 0.1)],
        color = ["red", "blue"],
        alpha = [0.1, 0.3]
    )
end

@testset "Tests from GeometryTypes" begin
    include("geometrytypes.jl")
end

@testset "Point & Vec type" begin
    include("fixed_arrays.jl")
end

end  # testset "GeometryBasics"


text:
@doc raw"""
```
albedo,illumination_direction,slant,tilt = estimate_img_properties(img::AbstractArray)
```
Attempts to calculate the illumination properties of an image for use in SFS
algorithms.
# Output
Returns albedo (ρ), illumination vector (I), slant (σ) and tilt (τ) of the image as Float,
Array{Float} of size 3, Float and Float respectively.
# Details
The illumination properties can be estimated using the first and second
moments and the gradients of the image (``E_x, E_y``) brightness as per the following:

```math
\begin{gathered}
\mu_1=\dfrac{\pi}{4}\rho\cos\sigma\\\\\mu_2=\dfrac{1}{6}\rho^2(1+3\cos^2\sigma)\\
\\\rho=\dfrac{\sqrt{6\pi\mu_2-48\mu_1^2}}{\pi}\\\\\cos\sigma=\dfrac{4\mu_1}{\sqrt
{6\pi\mu_2-48\mu_1^2}}\\\\\tan\tau=\dfrac{\text{average}(E_y)}{\text{average}
(E_x)}\\\\I=[\cos\tau\sin\sigma,\sin\tau\sin\sigma,\cos\sigma]
\end{gathered}
```
Note: Properties are an estimate and may not be correct under all conditions.
# Arguments
The function arguments are described in more detail below.
##  `img`
An `AbstractArray` storing the grayscale value of each pixel within
the range [0,1].
# Example
Compute properties of synthetic image generated using `generate_surface`.
```julia
using Images, ShapeFromShading

#generate synthetic image
img = generate_surface(SynthSphere(), 1, [0.2,0,0.9], radius = 5)

#estimate the properties
albedo, illumination_direction, slant, tilt = estimate_img_properties(img)
```
# Reference
1. T. Ping-Sing and M. Shah, "Shape from shading using linear approximation", Image and Vision Computing, vol. 12, no. 8, pp. 487-498, 1994. [doi:10.1016/0262-8856(94)90002-7](https://doi.org/10.1016/0262-8856(94)90002-7)
"""
function estimate_img_properties(img::AbstractArray)
    E = Array{Float64}(img)
    #calculate spatial gradient of img
    Ey, Ex = imgradients(E, KernelFactors.sobel, "replicate")
    #normalize gradient
    nEx = similar(E)
    nEy = similar(E)
    Exy = sqrt.((Ex.^2) .+ (Ey.^2))
    nEx = Ex ./ (Exy .+ eps())
    nEy = Ey ./ (Exy .+ eps())

    #calculate means
    μ₁ = mean(E)
    μ₂ = mean(E.^2)

    #calculate desired values
    g = sqrt((6.0 * (π^2) * μ₂) - (48.0 * (μ₁^2)))
    ρ = g / π
    t = (4 * μ₁) / g
    if t > 1.0
        t = 1.0
    elseif t < -1.0
        t = -1.0
    end
    σ = acos(t)
    τ = atan(mean(nEy) / mean(nEx))
    if τ < 0
        τ = τ + π;
    end
    I = [cos(τ) * sin(σ), sin(τ) * sin(σ), cos(σ)]
    return ρ, I, σ, τ
end


text:
@pyimport wntr.network.controls as controls
@pyimport wntr.sim.epanet as sims
function UpdateExtrema(wn::PyCall.PyObject, wn_dict::Dict{Any,Any}, Qmax_dict::Dict{String, Float64},
                        Qmin_dict::Dict{String,Float64}, Hmax_dict::Dict{String, Float64},
                        Hmin_dict::Dict{String,Float64}, num_sim::Int64, sim::Int64,
                        tight_coef::Float64, Q_lb::Float64)
    link_results, node_results = run_sims(wn)
    UpdateFlowRange(wn, wn_dict, link_results, node_results, Qmax_dict, Qmin_dict, num_sim, sim, tight_coef, Q_lb)
    UpdateHeadRange(wn_dict, Hmax_dict, Hmin_dict, link_results, node_results)
end

function run_sims(wn::PyCall.PyObject)
    wns = sims.EpanetSimulator(wn)
    results = wns[:run_sim]()
    link_results = results[:link]
    node_results = results[:node]
    return link_results, node_results
end
function UpdateFlowRange(wn::PyCall.PyObject, wn_dict::Dict{Any,Any},
                    link_results::Dict{Any,Any}, node_results::Dict{Any,Any},
                    Qmax_dict::Dict{String, Float64}, Qmin_dict::Dict{String,Float64},
                    num_sim::Int64, sim::Int64, tight_coef::Float64, Q_lb::Float64)

    pipes = wn_dict["pipes"]
    valves = wn_dict["valves"]
    links = merge(pipes, valves)
    for (key, link) in links
        name = link["name"]
        flows = link_results["flowrate"][name][:values]
        convert(Array{Float64}, flows)
        if sim == 1
            Qmin_dict[name] = convert(Float64,minimum(flows))
            Qmax_dict[name] = convert(Float64, maximum(flows))
        elseif sim < num_sim
            Qmin = Qmin_dict[name]
            Qmax = Qmax_dict[name]

            Qmin_dict[name] = minimum([Qmin, minimum(flows)])
            Qmax_dict[name] = maximum([Qmax, maximum(flows)])
        else
            Qmax_val = Qmax_dict[name]
            Qmin_val = Qmin_dict[name]
            gap = Qmax_val - Qmin_val
            if Qmin_val <= 0
                Qmin_dict[name] = Qmin_val - tight_coef*gap
            else
                Qmin_dict[name] = minimum([Q_lb, Qmin_val])
            end
            if Qmax_val >= 0
                Qmax_dict[name] = Qmax_val + tight_coef*gap
            else
                Qmax_dict[name] = maximum([Qmax_val, -Q_lb])
            end
        end
    end

    for (key, pump) in wn_dict["pumps"]
        name = pump["name"]
        flow = link_results["flowrate"][name][:values]
        convert(Array{Float64}, flow)
        if sim != 1
            flow = flow[flow .> 0]
        end

        if length(flow) !=0 || sim == num_sim
            if sim == 1
                Qmin_dict[name] = minimum(flow)
                Qmax_dict[name]= maximum(flow)
            elseif sim < num_sim
                Qmax = Qmax_dict[name]
                Qmin = Qmin_dict[name]
                Qmin_dict[name] = minimum([Qmin, minimum(flow)])
                Qmax_dict[name] = maximum([Qmax, maximum(flow)])
            else
                Qmin_val = Qmin_dict[name]
                Qmax_val = Qmax_dict[name]

                gap = Qmax_val - Qmin_val
                if pump["pump_type"] == "HEAD"
                    pump_curve_name = pump["pump_curve_name"]
                    pump_curve = wn_dict["curves"][pump_curve_name]["points"]
                else
                    pump_curve = [(pump["power"],0.0)] #power pump types gives fixed power value,
                end

                Qmin_dict[name] = maximum([pump_curve[1][1], Qmin_val - tight_coef*gap])
                Qmax_dict[name] = minimum([pump_curve[length(pump_curve)][1], Qmax_val + tight_coef*gap])
            end
        end
    end
end

function UpdateHeadRange(wn_dict::Dict{Any,Any}, Hmax_dict::Dict{String, Float64},
                        Hmin_dict::Dict{String,Float64}, link_results::Dict{Any,Any},
                        node_results::Dict{Any,Any})
    pumps = wn_dict["pumps"]
    nodes = wn_dict["nodes"]
    for (key, pump) in pumps
        name1 = pump["start_node_name"]
        elev1 = nodes[name1]["elevation"]

        name2 = pump["end_node_name"]
        elev2 = nodes[name2]["elevation"]

        pump_head = (node_results["pressure"][name2][:values] + elev2) - (node_results["pressure"][name1][:values] + elev1)
        pump_status = link_results["status"][key][:values]

        pump_head = pump_head .* pump_status #only care about head when pump is on
        pump_head = pump_head[pump_head .> 0]

        if length(pump_head) > 0
            try
                Hmin_dict[key] = min(minimum(pump_head), Hmin_dict[key])
                Hmax_dict[key] = max(maximum(pump_head), Hmax_dict[key])
            catch
                Hmin_dict[key] = minimum(pump_head)
                Hmax_dict[key] = maximum(pump_head)
            end
        end
    end
    return Hmin_dict, Hmax_dict
end


text:
#=
Interface for defining an evaluation policy
=#

"""
    evaluation(eval_policy, policy, env, obs, global_step, rng)
    returns the average reward of the current policy, the user can specify its own function
    f to carry the evaluation, we provide a default basic_evaluation that is just a rollout.
"""
function evaluation(f::Function, policy::AbstractNNPolicy, env::AbstractEnv, n_eval::Int64, max_episode_length::Int64, verbose::Bool = false)
    return f(policy, env, n_eval, max_episode_length, verbose)
end

# Examples
# just simulate the policy, return the average non discounted reward
function basic_evaluation(policy::AbstractNNPolicy, env::AbstractEnv, n_eval::Int64, max_episode_length::Int64, verbose::Bool)
    avg_r = 0.0
    avg_steps = 0.0
    for i=1:n_eval
        done = false
        r_tot = 0.0
        disc = 1.0
        step = 0
        reset!(env)
        obs = observe(env)
        ao = vcat(zeros(Float32, size((convert_a(Vector{Float32}, first(actions(env)), env.m)))), (obs))
        resetstate!(policy)
        while !done && step <= max_episode_length
            act = action(policy, ao)
            rew = act!(env, act)
            obs = observe(env)
            done = terminated(env)
            ao = vcat((convert_a(Vector{Float32}, act, env.m)), (obs))
            r_tot += disc * rew
            disc *= discount(env.m)
            step += 1
        end
        avg_steps += step
        avg_r += r_tot
    end
    if verbose
        @printf("Evaluation ... Avg Reward %2.2f | Avg Step %2.2f \n", avg_r/n_eval, avg_steps/n_eval)
    end
    return  avg_r / n_eval, avg_steps / n_eval, Dict()
end


text:
module FLOWVPM

import Dates
import GeometricTools: create_path

const RealFMM = Float64

# Available Kernels
const kernel_singular = (args...)->nothing
const kernel_gaussian = (args...)->nothing
const kernel_gaussianerf = (args...)->nothing
const kernel_winckelmans = (args...)->nothing
const singular = kernel_singular
const gaussian = kernel_gaussian
const gaussianerf = kernel_gaussianerf
const winckelmans = kernel_winckelmans

const rungekutta3 = (args...)->nothing
const euler = (args...)->nothing

mutable struct FMM
  # Optional user inputs
  p::Int32                        # Multipole expansion order
  ncrit::Int32                    # Max number of particles per leaf
  theta::RealFMM                  # Neighborhood criterion
  phi::RealFMM                    # Regularizing neighborhood criterion

  FMM(; p=4, ncrit=50, theta=0.4, phi=1/3) = new(p, ncrit, theta, phi)
end

"""
    `ViscousScheme{R}`

Type declaring viscous scheme.

Implementations must have the following properties:
    * `nu::R`                   : Kinematic viscosity.
"""
abstract type ViscousScheme{R} end

"""
Implementation of viscous diffusion scheme that gets called in the inner loop
of the time integration scheme at each time step.
"""
function viscousdiffusion(pfield, scheme::ViscousScheme, dt; optargs...)
    error("Viscous diffusion scheme has not been implemented yet!")
end

viscousdiffusion(pfield, dt; optargs...
                    ) = viscousdiffusion(pfield, pfield.viscous, dt; optargs...)

zeta_direct(args...) = nothing
zeta_fmm(args...) = nothing
##### END OF ABSTRACT VISCOUS SCHEME ###########################################

################################################################################
# INVISCID SCHEME TYPE
################################################################################
struct Inviscid{R} <: ViscousScheme{R}
    nu::R                                 # Kinematic viscosity
    Inviscid{R}(; nu=zero(R)) where {R} = new(nu)
end

Inviscid() = Inviscid{RealFMM}()

"""
    `isinviscid(scheme::ViscousScheme)`

Returns true if viscous scheme is inviscid.
"""
isinviscid(scheme::ViscousScheme) = typeof(scheme).name == Inviscid.body.name

viscousdiffusion(pfield, scheme::Inviscid, dt; optargs...) = nothing
##### END OF INVISCID SCHEME ###################################################

################################################################################
# CORE SPEADING SCHEME TYPE
################################################################################
mutable struct CoreSpreading{R} <: ViscousScheme{R}
    # User inputs
    nu::R                                 # Kinematic viscosity
    sgm0::R                               # Core size after reset
    zeta::Function                        # Basis function evaluation method

    # Optional inputs
    beta::R                               # Maximum core size growth σ/σ_0
    itmax::Int                            # Maximum number of RBF iterations
    tol::R                                # RBF interpolation tolerance
    iterror::Bool                         # Throw error if RBF didn't converge
    verbose::Bool                         # Verbose on RBF interpolation
    v_lvl::Int                            # Verbose printing tab level
    debug::Bool                           # Print verbose for debugging

    # Internal properties
    t_sgm::R                              # Time since last core size reset
    rbf::Function                         # RBF function
    rr0s::Array{R, 1}                     # Initial field residuals
    rrs::Array{R, 1}                      # Current field residuals
    prev_rrs::Array{R, 1}                 # Previous field residuals
    pAps::Array{R, 1}                     # pAp product
    alphas::Array{R, 1}                   # Alpha coefficients
    betas::Array{R, 1}                    # Beta coefficients
    flags::Array{Bool, 1}                 # Convergence flags

    CoreSpreading{R}(
                        nu, sgm0, zeta;
                        beta=R(1.5),
                        itmax=R(15), tol=R(1e-3),
                        iterror=true, verbose=false, v_lvl=2, debug=false,
                        t_sgm=R(0.0),
                        rbf=rbf_conjugategradient,
                        rr0s=zeros(R, 3), rrs=zeros(R, 3), prev_rrs=zeros(R, 3),
                        pAps=zeros(R, 3), alphas=zeros(R, 3), betas=zeros(R, 3),
                        flags=zeros(Bool, 3)
                    ) where {R} = new(
                        nu, sgm0, zeta,
                        beta,
                        itmax, tol,
                        iterror, verbose, v_lvl, debug,
                        t_sgm,
                        rbf,
                        rr0s, rrs, prev_rrs,
                        pAps, alphas, betas,
                        flags
                    )
end

CoreSpreading(nu, sgm0, args...; optargs...
                    ) = CoreSpreading{RealFMM}(RealFMM(nu), RealFMM(sgm0), args...; optargs...)

"""
   `iscorespreading(scheme::ViscousScheme)`

Returns true if viscous scheme is core spreading.
"""
iscorespreading(scheme::ViscousScheme
                            ) = typeof(scheme).name == CoreSpreading.body.name

mutable struct ParticleField{T, V<:ViscousScheme}
    # User inputs
    maxparticles::Int                           # Maximum number of particles
    particles                                   # Array of particles
    bodies                                      # ExaFMM array of bodies
    viscous::V                                  # Viscous scheme

    # Internal properties
    np::Int                                     # Number of particles in the field
    nt::Int                                     # Current time step number
    t::Float64                                  # Current time

    # Solver setting
    kernel                                      # Vortex particle kernel
    UJ::Function                                # Particle-to-particle calculation

    # Optional inputs
    Uinf::Function                              # Uniform freestream function Uinf(t)
    transposed::Bool                            # Transposed vortex stretch scheme
    relax::Bool                                 # Activates relaxation scheme
    rlxf::Float64                               # Relaxation factor (fraction of dt)
    integration::Function                       # Time integration scheme
    fmm::FMM                                    # Fast-multipole settings

    ParticleField{T, V}(
                        maxparticles,
                        particles, bodies, viscous;
                        np=0, nt=0, t=0.0,
                        kernel=gaussianerf,
                        UJ=UJ_fmm,
                        Uinf=t->zeros(3),
                        transposed=true,
                        relax=true, rlxf=0.3,
                        integration=rungekutta3,
                        fmm=FMM(),
                 ) where {T, V} = new(
                        maxparticles,
                        particles, bodies, viscous,
                        np, nt, t,
                        kernel,
                        UJ,
                        Uinf,
                        transposed,
                        relax, rlxf,
                        integration,
                        fmm,
                  )
end

function ParticleField(maxparticles; viscous::V=Inviscid(),
                                            optargs...) where {V<:ViscousScheme}
    return ParticleField{RealFMM, Inviscid{RealFMM}}(maxparticles, nothing, nothing, viscous; optargs...)
end

function get_np(self::ParticleField)
  return self.np
end
function add_particle(args...; optargs...)
    nothing
end
function remove_particle(args...; optargs...)
    nothing
end
function UJ_direct(args...; optargs...)
    nothing
end
function UJ_fmm(args...; optargs...)
    nothing
end
function get_particleiterator(args...; optargs...)
    nothing
end
iterator(args...; optargs...) = get_particleiterator(args...; optargs...)
iterate(args...; optargs...) = get_particleiterator(args...; optargs...)
get_X(self::ParticleField, i::Int) = nothing
get_Gamma(self::ParticleField, i::Int) = nothing
get_sigma(self::ParticleField, i::Int) = nothing
get_U(self::ParticleField, i::Int) = nothing
get_W(self::ParticleField, i::Int) = nothing
isinviscid(self::ParticleField) = isinviscid(self.viscous)

"""
  `run_vpm!(pfield, dt, nsteps; runtime_function=nothing, save_path=nothing,
run_name="pfield", nsteps_save=1, verbose=true, prompt=true)`

Solves `nsteps` of the particle field with a time step of `dt`.

**Optional Arguments**
* `runtime_function::Function`   : Give it a function of the form
                            `myfun(pfield, t, dt)`. On each time step it
                            will call this function. Use this for adding
                            particles, deleting particles, etc.
* `static_particles_function::Function`   : Give it a function of the form
                            `myfun(pfield, t, dt)` to add static particles
                            representing solid boundaries to the solver. This
                            function is called at every time step right before
                            solving the governing equations, and any new
                            particles added by this function are immediately
                            removed.
* `nsteps_relax::Int`   : Relaxes the particle field every this many time steps.
* `save_path::String`   : Give it a string for saving VTKs of the particle
                            field. Creates the given path.
* `run_name::String`    : Name of output files.
* `nsteps_save::Int64`  : Saves vtks every this many time steps.
* `prompt::Bool`        : If `save_path` already exist, it will prompt the
                            user before overwritting the folder if true; it will
                            directly overwrite it if false.
* `verbose::Bool`       : Prints progress of the run to the terminal.
* `verbose_nsteps::Bool`: Number of time steps between verbose.
"""
function run_vpm!(pfield::ParticleField, dt::Real, nsteps::Int;
                      # RUNTIME OPTIONS
                      runtime_function::Function=(pfield, t, dt)->false,
                      # OUTPUT OPTIONS
                      save_path::Union{Nothing, String}=nothing,
                      create_savepath::Bool=true,
                      run_name::String="pfield",
                      save_code::String="",
                      prompt::Bool=true,
                      verbose::Bool=true, verbose_nsteps::Int=10, v_lvl::Int=0,
                      optargs...)

    # Creates save path and save code
    if save_path!=nothing && create_savepath
        create_path(save_path, prompt)
    end
    # Save code
    if save_path!=nothing && save_code!=""
        cp(save_code, save_path*"/"; force=true)
    end

    run_id = save_path!=nothing ? joinpath(save_path, run_name) : ""

    if verbose
        time_beg = Dates.DateTime(Dates.now())
        println("\t"^v_lvl*"*"^(73-8*v_lvl)*"\n"*"\t"^v_lvl*"START $run_id\t$time_beg\n"*
                "\t"^v_lvl*"*"^(73-8*v_lvl))
    end

    # RUN
    for i in 0:nsteps

        if verbose && i%verbose_nsteps==0
            println("\t"^(v_lvl+1)*"Time step $i out of $nsteps"*
            "\tParticles: $(get_np(pfield))")
        end

        # Time step
        if i!=0
            pfield.nt +=1
            pfield.t += dt
        end

        # Calls user-defined runtime function
        breakflag = runtime_function(pfield, pfield.t, dt)

        # User-indicated end of simulation
        if breakflag
            break
        end

    end

    if verbose
        time_end = Dates.DateTime(Dates.now())
        hrs,mins,secs = timeformat(time_beg, time_end)
        println("\t"^v_lvl*"*"^(73-8*v_lvl))
        println("\t"^v_lvl*"END $run_id\t$time_end")
        println("\t"^v_lvl*"*"^(73-8*v_lvl))
        println("\t"^v_lvl*"ELAPSED TIME: $hrs hours $mins minutes $secs seconds")
    end

    return nothing
end

function timeformat(time_beg, time_end)
    time_delta = Dates.value(time_end)-Dates.value(time_beg)
    hrs = Int(floor(time_delta/1000/60/60))
    mins = Int(floor((time_delta-hrs*60*60*1000)/1000/60))
    secs = Int(floor((time_delta-hrs*60*60*1000-mins*1000*60)/1000))
    return hrs,mins,secs
end

end # END OF MODULE

max_stars_repo_name: dendisuhubdy/FLOWUnsteady

max_issues_repo_name: dendisuhubdy/FLOWUnsteady

max_forks_repo_name: dendisuhubdy/FLOWUnsteady


text:
# ---
# title: 1041. Robot Bounded In Circle
# id: problem1041
# author: Tian Jun
# date: 2020-10-31
# difficulty: Medium
# categories: Math
# link: <https://leetcode.com/problems/robot-bounded-in-circle/description/>
# hidden: true
# ---
#
# On an infinite plane, a robot initially stands at `(0, 0)` and faces north.
# The robot can receive one of three instructions:
#
#   * `"G"`: go straight 1 unit;
#   * `"L"`: turn 90 degrees to the left;
#   * `"R"`: turn 90 degress to the right.
#
# The robot performs the `instructions` given in order, and repeats them
# forever.
#
# Return `true` if and only if there exists a circle in the plane such that the
# robot never leaves the circle.
#
#
#
# **Example 1:**
#
#
#
#     Input: "GGLLGG"
#     Output: true
#     Explanation:
#     The robot moves from (0,0) to (0,2), turns 180 degrees, and then returns to (0,0).
#     When repeating these instructions, the robot remains in the circle of radius 2 centered at the origin.
#
#
# **Example 2:**
#
#
#
#     Input: "GG"
#     Output: false
#     Explanation:
#     The robot moves north indefinitely.
#
#
# **Example 3:**
#
#
#
#     Input: "GL"
#     Output: true
#     Explanation:
#     The robot moves from (0, 0) -> (0, 1) -> (-1, 1) -> (-1, 0) -> (0, 0) -> ...
#
#
#
#
# **Note:**
#
#   1. `1 <= instructions.length <= 100`
#   2. `instructions[i]` is in `{'G', 'L', 'R'}`
#
#
## @lc code=start
using LeetCode

## add your code here:
## @lc code=end


text:
#############
# Repo Type #
#############

type Repo <: GitLabType
    name::Nullable{GitLabString}
    visibility_level::Nullable{Int}
    homepage::Nullable{HttpCommon.URI}
    git_http_url::Nullable{HttpCommon.URI}
    description::Nullable{GitLabString}
    project_id::Nullable{Int}

    id::Nullable{Int}
    default_branch::Nullable{GitLabString}
    tag_list::Nullable{Vector{GitLabString}}
    public::Nullable{Bool}
    archived::Nullable{Bool}
    ## TODO FIX ssh_url_to_repo::Nullable{HttpCommon.URI}
    http_url_to_repo::Nullable{HttpCommon.URI}
    web_url::Nullable{HttpCommon.URI}
    owner::Nullable{Owner}
    name_with_namespace::Nullable{GitLabString}
    path::Nullable{GitLabString}
    path_with_namespace::Nullable{GitLabString}
    issues_enabled::Nullable{Bool}
    merge_requests_enabled::Nullable{Bool}
    wiki_enabled::Nullable{Bool}
    builds_enabled::Nullable{Bool}
    snippets_enabled::Nullable{Bool}
    container_registry_enabled::Nullable{Bool}
    created_at::Nullable{Dates.DateTime}
    last_activity_at::Nullable{Dates.DateTime}
    shared_runners_enabled::Nullable{Bool}
    creator_id::Nullable{Int}
    ## TODO FIX namespace::Nullable{Namespace}
    ## \"namespace\":{\"id\":2,\"name\":\"mdpradeep\",\"path\":\"mdpradeep\",\"owner_id\":2,\"created_at\":\"2016-06-17T07:09:56.494Z\",\"updated_at\":\"2016-06-17T07:09:56.494Z\",\"description\":\"\",\"avatar\":null,\"share_with_group_lock\":false,\"visibility_level\":20}
    avatar_url::Nullable{HttpCommon.URI}
    star_count::Nullable{Int}
    forks_count::Nullable{Int}
    open_issues_count::Nullable{Int}
    runners_token::Nullable{GitLabString}
    public_builds::Nullable{Bool}
    ## TODO permissions::Nullable{Permissions}
    ## \"permissions\":{\"project_access\":{\"access_level\":40,\"notification_level\":3},\"group_access\":null}
end

Repo(data::Dict) = json2gitlab(Repo, data)
## MDP Repo(full_name::AbstractString) = Repo(Dict("full_name" => full_name))
Repo(full_name::AbstractString) = Repo(Dict("name" => full_name))

## MDP namefield(repo::Repo) = repo.full_name
namefield(repo::Repo) = repo.name

###############
# API Methods #
###############

# repos #
#-------#

function repo_by_name(repo_name; options...)
    result = gh_get_json("/api/v3/projects/search/$(repo_name)"; options...)
    return Repo(result[1])
end

function repo(id; options...)
    ## result = gh_get_json("/repos/$(name(repo_obj))"; options...)
    result = gh_get_json("/api/v3/projects/$(id)"; options...)
    return Repo(result)
end

# forks #
#-------#

function forks(repo; options...)
    error("Not implemented yet !!")
    ## TODO
    ## results, page_data = gh_get_paged_json("/repos/$(name(repo))/forks"; options...)
    results, page_data = gh_get_paged_json("/api/v3/projects/$(get(repo.id))/forks"; options...)
    return map(Repo, results), page_data
end

function create_fork(repo; options...)
    ## result = gh_post_json("/repos/$(name(repo))/forks"; options...)
    result = gh_post_json("/api/v3/projects/fork/$(get(repo.id))"; options...)
    return Repo(result)
end

function delete_fork(repo; options...)
    ## /projects/:id/fork
    result = gh_delete_json("/api/v3/projects/$(get(repo.id))/fork"; options...)
    return Repo(result)
end

# contributors/collaborators #
#----------------------------#

function contributors(repo; options...)
    ## results, page_data = gh_get_paged_json("/repos/$(name(repo))/contributors"; options...)
    results, page_data = gh_get_paged_json("/api/v3/projects/$(get(repo.id))/repository/contributors"; options...)
    results = [Dict("contributor" => Owner(i), "contributions" => i["commits"]) for i in results]
    return results, page_data
end

function collaborators(repo; options...)
    ## MDP results, page_data = gh_get_json("/repos/$(name(repo))/collaborators"; options...)
    ## MDP results, page_data = gh_get_json("/api/v3/projects/$(get(repo.id))/repository/contributors"; options...)
    ## results, page_data = gh_get_paged_json("/api/v3/projects/$(get(repo.id))/members"; options...)
    results = gh_get_json("/api/v3/projects/$(get(repo.id))/members"; options...)
    return map(Owner, results)
end

function iscollaborator(repo, user; options...)
    collaborators = GitLab.collaborators(repo; options...)
    for c in collaborators
        if get(c.username) == user
            return true
        end
    end

    return false
end

function add_collaborator(repo, user; options...)
    ## MDP path = "/repos/$(name(repo))/collaborators/$(name(user))"
    ## path = "/api/v3/projects/$(get(repo.id))/members/$(user)"
    ## return gh_put(path; options...)
    path = "/api/v3/projects/$(get(repo.id))/members"
    return gh_post(path; options...)
end

function remove_collaborator(repo, user; options...)
    ## MDP path = "/repos/$(name(repo))/collaborators/$(name(user))"
    ## path = "/api/v3/projects/$(get(repo.id))/repository/contributors/$(name(user))"
    path = "/api/v3/projects/$(get(repo.id))/members/$(user)"
    return gh_delete(path; options...)
end

# stats #
#-------#

## TODO Check how to enable sidekiq stats !
function stats(repo, stat, attempts = 3; options...)
    ## MDP path = "/repos/$(name(repo))/stats/$(name(stat))"
    ## path = "/api/v3/projects/$(get(repo.id))/repository/stats/$(name(stat))"
    path = "/api/v3/projects/sidekiq/$(name(stat))"
    local r
    for a in 1:attempts
        r = gh_get(path; handle_error = false, options...)
        r.status == 200 && return r
        sleep(2.0)
    end
    return r
end


text:
# ---
# title: 1051. Height Checker
# id: problem1051
# author: Tian Jun
# date: 2020-10-31
# difficulty: Easy
# categories: Array
# link: <https://leetcode.com/problems/height-checker/description/>
# hidden: true
# ---
#
# Students are asked to stand in non-decreasing order of heights for an annual
# photo.
#
# Return the minimum number of students that must move in order for all students
# to be standing in non-decreasing order of height.
#
# Notice that when a group of students is selected they can reorder in any
# possible way between themselves and the non selected students remain on their
# seats.
#
#
#
# **Example 1:**
#
#
#
#     Input: heights = [1,1,4,2,1,3]
#     Output: 3
#     Explanation:
#     Current array : [1,1,4,2,1,3]
#     Target array  : [1,1,1,2,3,4]
#     On index 2 (0-based) we have 4 vs 1 so we have to move this student.
#     On index 4 (0-based) we have 1 vs 3 so we have to move this student.
#     On index 5 (0-based) we have 3 vs 4 so we have to move this student.
#
#
# **Example 2:**
#
#
#
#     Input: heights = [5,1,2,3,4]
#     Output: 5
#
#
# **Example 3:**
#
#
#
#     Input: heights = [1,2,3,4,5]
#     Output: 0
#
#
#
#
# **Constraints:**
#
#   * `1 <= heights.length <= 100`
#   * `1 <= heights[i] <= 100`
#
#
## @lc code=start
using LeetCode

## add your code here:
## @lc code=end


text:
# *****************************************************************************
# Written by Ritchie Lee, ritchie.lee@sv.cmu.edu
# *****************************************************************************
# Copyright ã 2015, United States Government, as represented by the
# Administrator of the National Aeronautics and Space Administration. All
# rights reserved.  The Reinforcement Learning Encounter Simulator (RLES)
# platform is licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License. You
# may obtain a copy of the License at
# http://www.apache.org/licenses/LICENSE-2.0. Unless required by applicable
# law or agreed to in writing, software distributed under the License is
# distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
# KIND, either express or implied. See the License for the specific language
# governing permissions and limitations under the License.
# _____________________________________________________________________________
# Reinforcement Learning Encounter Simulator (RLES) includes the following
# third party software. The SISLES.jl package is licensed under the MIT Expat
# License: Copyright (c) 2014: Youngjun Kim.
# Permission is hereby granted, free of charge, to any person obtaining a copy
# of this software and associated documentation files (the "Software"), to
# deal in the Software without restriction, including without limitation the
# rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
# sell copies of the Software, and to permit persons to whom the Software is
# furnished to do so, subject to the following conditions:
# The above copyright notice and this permission notice shall be included in
# all copies or substantial portions of the Software. THE SOFTWARE IS PROVIDED
# "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT
# NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR
# PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
# HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN
# ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
# CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
# *****************************************************************************

"""
Generic decision tree based on callbacks
"""
module DecisionTrees

export DTParams, DecisionTree, DTNode, build_tree, classify, isleaf, get_max_depth,
    get_members, count_leafs, count_nonleafs
export get_truth, get_splitter, get_split_labels #to be overridden

using StatsBase
using RLESUtils, TreeUtils, TreeIterators
import RLESTypes: SymbolTable
import Base: length
import TreeIterators.get_children

type DTParams
    num_data::Int64 #number of datapoints in training set
    maxdepth::Int64
    predict_type::Type #T1
    label_type::Type #T2
    userargs::SymbolTable #additional variables in callabacks
    id::Int64 #hack, use this to store state of node id's assigned, init to 0

    DTParams(num_data::Int64, maxdepth::Int64, predict_type::Type, label_type::Type,
        userargs::SymbolTable) = new(num_data, maxdepth, predict_type, label_type,
        userargs, 0) #initialize id=0
end

type DTNode{T1,T2}
    id::Int64
    depth::Int64
    members::Vector{Int64} #indices into data starting at 1
    split_rule::Any #object used in callback for split rule
    children::Dict{T1,DTNode{T1,T2}} #key=split_rule predicts, value=child node, T1=predict label type
    label::T2 #T2=label type
    confidence::Float64
end

type DecisionTree{T1,T2}
    root::DTNode{T1,T2}
end

#User should override these functions
get_truth() = error("get_truth not defined")
get_splitter() = error("get_splitter not defined")
get_split_labels() = error("get_split_labels not defined")

get_children(node::DTNode) = collect(values(node.children))

function build_tree(p::DTParams)
    members = collect(1:p.num_data)
    depth = 0
    root = process_child(p, members, depth)
    DecisionTree(root)
end

function process_child(p::DTParams, members::Vector{Int64}, depth::Int64)
    p.id += 1 #node id's are ordered depth-first
    id = p.id
    members_copy = deepcopy(members)
    labels = get_truth(members, p.userargs)
    label = mode(labels)
    confidence = proportions(labels, label:label)[1] #proportion of times label appears
    T1 = p.predict_type
    T2 = p.label_type
    @assert issubtype(eltype(labels), T2)
    children = Dict{T1,DTNode{T1,T2}}()
    if confidence == 1.0 || depth >= p.maxdepth #all outputs same, or reached stopping criterion
        split_rule = nothing
    else #consider splitting
        split_rule = get_splitter(members, p.userargs)
        if split_rule != nothing #get_splitter returns nothing if can't find a good split
            predict_labels = get_split_labels(split_rule, members, p.userargs)
            @assert issubtype(eltype(predict_labels), T1)
            predict_set = sort(unique(predict_labels))
            child_nodes = map(predict_set) do label_class
                #TODO: parallelize, pmap crashes due to missing callbacks
                child_members = members[find(x -> x == label_class, predict_labels)]
                child = process_child(p, child_members, depth + 1)
                return child
            end
            for (label_class, child_node) in zip(predict_set, child_nodes)
                children[label_class] = child_node
            end
        end
    end
    DTNode{T1,T2}(id, depth, members_copy, split_rule, children, label, confidence)
end

classify{T1,T2}(p::DTParams, tree::DecisionTree{T1,T2}, x) =
    classify(p, tree.root, x) #x could be either id or data

function classify{T1,T2}(p::DTParams, node::DTNode{T1,T2}, x)
    if isempty(node.children) #leaf node
        return node.label::T2
    end
    label = get_split_labels(node.split_rule, [x], p.userargs)[1]
    child = node.children[label]
    return classify(p, child, x)
end

isleaf{T1,T2}(node::DTNode{T1,T2}) = isempty(node.children)

get_max_depth{T1,T2}(dtree::DecisionTree{T1,T2}) = get_max_depth(dtree.root)
function get_max_depth{T1,T2}(node::DTNode{T1,T2})
    if isleaf(node)
        return node.depth
    else
        return map(get_max_depth, values(node.children)) |> maximum
    end
end

get_members{T1,T2}(tree::DecisionTree{T1,T2}) = get_members(tree.root)
get_members{T1,T2}(node::DTNode{T1,T2}) = node.members

count_nonleafs{T1,T2}(tree::DecisionTree{T1,T2}) =
    count_nodes(x->!DecisionTrees.isleaf(x), tree.root)

count_leafs{T1,T2}(tree::DecisionTree{T1,T2}) =
    count_nodes(DecisionTrees.isleaf, tree.root)

length{T1,T2}(tree::DecisionTree{T1,T2}) = count_nodes(tree.root)

end #module