text:
##################################################################################
# This file is part of StateSpaceEcon.jl
# BSD 3-Clause License
# Copyright (c) 2020-2022, Bank of Canada
# All rights reserved.
##################################################################################

"""
    Plans

Module part of StateSpaceEcon. This module implements the `Plan` data structure,
which is used in simulations. The plan object contains information about the
range of the simulation and which variables and shocks are exogenous or
endogenous at each period of the range.

### Constructors
  * [`Plan(model, range)`](@ref Plan)

### Modify the plan
  * [`exogenize!`](@ref), [`endogenize!`](@ref) - make variables exogenous or
    endogenous
  * [`exog_endo!`](@ref), [`endo_exog!`](@ref) - swap exogenous and endogenous
    variables
  * [`autoexogenize!`](@ref) - exogenize and endogenize variables according to
    the list in the model

### Prepare data for simulation
  * [`zeroarray`](@ref), [`zerodict`](@ref), [`zerodata`](@ref) - prepare a
    matrix or a dictionary or a [`SimData`](@ref) of data for the simulation
    containing zeros.
  * [`steadystatearray`](@ref), [`steadystatedict`](@ref),
    [`steadystatedata`](@ref) - prepare a matrix or a dictionary or a [`SimData`](@ref) of data for the
    simulation containing the steady state.

"""
module Plans

using TimeSeriesEcon
using ModelBaseEcon

# exported from this module that are also exported from StateSpaceEcon
export Plan,
    exogenize!, endogenize!, exog_endo!, endo_exog!, autoexogenize!,
    zeroarray, zerodict, steadystatearray, steadystatedict, zerodata, steadystatedata

# exported from this module, but meant for internal use in StateSpaceEcon
export plansum

"""
    Plan{T <: MIT}

A data structure representing the simulation plan. It holds information about
the time range of the simulation and which variables/shocks are exogenous at
each period.

"""
struct Plan{T<:MIT} <: AbstractVector{Vector{Symbol}}
    range::AbstractUnitRange{T}
    varshks::NamedTuple
    exogenous::BitArray{2}
end

Plan(model::Model, r::Union{MIT,Int}) = Plan(model, r:r)

"""
    Plan(model, range)

Create a default simulation plan for the given model over the given range. The
range of the plan is augmented to include periods before and after the given
range, over which initial and final conditions will be applied.

Instead of a range, one could also pass in a single moment in time
([`MIT`](@ref TimeSeriesEcon.MIT)) instance, in which case it is interpreted as
a range of length 1.

"""
function Plan(model::Model, range::AbstractUnitRange)
    if !(eltype(range) <: MIT)
        range = UnitRange{MIT{Unit}}(range)
    end
    range = (first(range)-model.maxlag):(last(range)+model.maxlead)
    local varshks = model.varshks
    local N = length(varshks)
    local names = tuple(Symbol.(varshks)...) # force conversion of Vector{ModelSymbol} to NTuple{N,Symbol}
    p = Plan(range, NamedTuple{names}(1:N), falses(length(range), length(varshks)))
    for (ind, var) = enumerate(varshks)
        if isshock(var) || isexog(var)
            p.exogenous[:, ind] .= [true]
        end
    end
    return p
end

#######################################
# AbstractVector interface

@inline Base.size(p::Plan) = size(p.range)
@inline Base.axes(p::Plan) = (p.range,)
@inline Base.length(p::Plan) = length(p.range)
@inline Base.IndexStyle(::Plan) = IndexLinear()
@inline Base.similar(p::Plan) = Plan(p.range, p.varshks, similar(p.exogenous))
@inline Base.copy(p::Plan) = Plan(p.range, p.varshks, copy(p.exogenous))

@inline _offset(p::Plan{T}, idx::T) where {T<:MIT} = convert(Int, idx - first(p.range) + 1)
@inline _offset(p::Plan{T}, idx::AbstractUnitRange{T}) where {T<:MIT} =
    _offset(p, first(idx)):_offset(p, last(idx))

## Integer index is taken as is
@inline Base.getindex(p::Plan, idx::Int) = p[p.range[idx]]
@inline Base.getindex(p::Plan, idx::AbstractUnitRange{Int}) = p[p.range[idx]]
@inline Base.getindex(p::Plan, idx::Int, ::Val{:inds}) = findall(p.exogenous[idx, :])

# Index of the frequency type returns the list of exogenous symbols
@inline Base.getindex(p::Plan{T}, idx::T) where {T<:MIT} = [keys(p.varshks)[p[_offset(p, idx), Val(:inds)]]...,]
@inline function Base.getindex(p::Plan{T}, idx::T, ::Val{:inds}) where {T<:MIT}
    first(p.range) ≤ idx ≤ last(p.range) || throw(BoundsError(p, idx))
    return p[_offset(p, idx), Val(true)]
end

# A range returns the plan trimmed over that exact range.
@inline Base.getindex(p::Plan{MIT{Unit}}, rng::AbstractUnitRange{Int}) = p[UnitRange{MIT{Unit}}(rng)]
@inline function Base.getindex(p::Plan{T}, rng::AbstractUnitRange{T}) where {T<:MIT}
    rng.start < p.range.start && throw(BoundsError(p, rng.start))
    rng.stop > p.range.stop && throw(BoundsError(p, rng.stop))
    return Plan{T}(rng, p.varshks, p.exogenous[_offset(p, rng), :])
end

# A range with a model returns a plan trimmed over that range and extended for initial and final conditions.
@inline Base.getindex(p::Plan{MIT{Unit}}, rng::AbstractUnitRange{Int}, m::Model) = p[UnitRange{MIT{Unit}}(rng), m]
@inline function Base.getindex(p::Plan{T}, rng::AbstractUnitRange{T}, m::Model) where {T<:MIT}
    rng = (rng.start-m.maxlag):(rng.stop+m.maxlead)
    return p[rng]
end

@inline Base.setindex!(p::Plan, x, i...) = error("Cannot assign directly. Use `exogenize` and `endogenize` to alter plan.")

#######################################
# query the exo-end status of a variable

@inline Base.getindex(p::Plan{T}, vars::Symbol...) where {T} = begin
    var_inds = [p.varshks[v] for v in vars]
    Plan{T}(p.range, NamedTuple{(vars...,)}(eachindex(vars)), p.exogenous[:, var_inds])
end

#######################################
# Pretty printing

@inline Base.summary(io::IO, p::Plan) = print(io, typeof(p), " with range ", p.range)

# Used in the show() implementation below
function collapsed_range(p::Plan{T}) where {T<:MIT}
    ret = Pair{Union{T,UnitRange{T}},Vector{Symbol}}[]
    i1 = first(p.range)
    i2 = i1
    val = p[i1]
    @inline make_key() = i1 == i2 ? i1 : i1:i2
    for i in p.range[2:end]
        val_i = p[i]
        if val_i == val
            i2 = i
        else
            push!(ret, make_key() => val)
            i1 = i2 = i
            val = val_i
        end
    end
    push!(ret, make_key() => val)
end

@inline Base.show(io::IO, ::MIME"text/plain", p::Plan) = show(io, p)
function Base.show(io::IO, p::Plan)
    # 0) show summary before setting :compact
    summary(io, p)
    isempty(p) && return
    # print(io, ":")
    nrow, ncol = displaysize(io)
    limit = get(io, :limit, true)
    cp = collapsed_range(p)
    # find the longest string left of "=>" for padding
    maxl = maximum(length("$k") for (k, v) in cp)
    if limit
        dcol = ncol - maxl - 6
    else
        dcol = typemax(Int)
    end
    function print_exog(names)
        if isempty(names)
            print(io, "∅")
            return
        end
        lens = cumsum(map(x -> length("$x, "), names))
        show = lens .< dcol
        show[1] = true
        print(io, join(names[show], ", "))
        if !all(show)
            print(io, ", …")
        end
    end
    if !limit || length(cp) <= nrow - 5
        for (r, v) in cp
            print(io, "\n  ", lpad("$r", maxl, " "), " → ")
            print_exog(v)
        end
    else
        top = div(nrow - 5, 2)
        bot = length(cp) - nrow + 6 + top
        for (r, v) in cp[1:top]
            print(io, "\n  ", lpad("$r", maxl, " "), " → ")
            print_exog(v)
        end
        print(io, "\n   ⋮")
        for (r, v) in cp[bot:end]
            print(io, "\n  ", lpad("$r", maxl, " "), " → ")
            print_exog(v)
        end
    end
end

#######################################
# export and import Plan instances

export exportplan, importplan

@inline exportplan(p::Plan; kwargs...) = exportplan(Base.stdout, p; kwargs...)
@inline exportplan(file::AbstractString, p::Plan; kwargs...) = (
    open(file, "w") do f
        exportplan(f, p; kwargs...)
    end
)

# center padding
_cpad(x, n::Integer, args...) = _cpad(string(x), Int(n), args...)
function _cpad(x::AbstractString, n::Integer, args...)
    lx = length(x)
    return signed(n) < 2 + lx ? lpad(x, n, args...) :
           lpad(rpad(x, (signed(n) + lx) ÷ 2, args...), n, args...)
end

function exportplan(io::IO, p::Plan;
    alphabetical = false,  # whether to sort variables
    exog_mark = "X", endo_mark = "-",  # symbols used for each class
    delim = " ",  # set to "," to get a CSV file (with 3 skip rows and 1 header row)
    _name_delim = delim,  # padding after NAME column
    _range_delim = delim    # padding between range columns
)
    summary(io, p)
    println(io)
    println(io, "Range: ", p.range)
    println(io, "Variables: ", p.varshks)
    width1 = 2 + maximum(length, (sprint(print, v; context = io, sizehint = 20) for v in keys(p.varshks)))
    width1 = max(width1, 2 + length("NAME"))
    ranges, _ = zip(collapsed_range(p)...)
    width2 = 1 .+ map(length, (sprint(print, rng; context = io, sizehint = 15) for rng in ranges))
    width2 = max.(width2, 1 + maximum(length, (exog_mark, endo_mark)))
    tf_matr = p.exogenous[map(x -> _offset(p, first(x)), [ranges...]), :]
    println(io, "($(exog_mark)) = Exogenous, ($(endo_mark)) = Endogenous:")
    print(io, lpad("NAME", width1), _name_delim)
    tmp = (lpad(rng, w) for (rng, w) in zip(ranges, width2))
    println(io, join(tmp, _range_delim))
    varind = pairs(p.varshks)
    if alphabetical
        varind = sort([varind...], by = Base.Fix2(getindex, 1),
            lt = (l, r) -> isless(string(l), string(r)))
    end
    for (var, ind) in varind
        print(io, lpad(var, width1), _name_delim)
        tmp = (_cpad(b ? exog_mark : endo_mark, w) for (w, b) in zip(width2, tf_matr[:, ind]))
        println(io, join(tmp, _range_delim))
    end
end

@inline function importplan(fname::AbstractString; kwargs...)
    open(fname, "r") do f
        importplan(f; kwargs...)
    end
end

function importplan(io::IO)
    # parse line 1. Example: "Plan{MIT{Quarterly}} with range 2000Q1:2100Q1"
    line = readline(io)
    m = match(r"Plan\{MIT\{(\w+)\}\} with range ([\w:]+)", line)
    if m === nothing
        error("expected Plan{MIT{Frequency}} at the start of line 1, got ", line, ".")
    end
    freq = parse_frequency(m.captures[1], 1)
    rng = parse_range(m.captures[2], false, 1)
    # parse line 2. Example: "Range: 2000Q1:2100Q1"  range must match line 1
    line = readline(io)
    m = match(r"Range: ([\w:]+)", line)
    if m === nothing
        error("expected Range: at the start of line 2, got ", line, ".")
    end
    rng1 = parse_range(m.captures[1], false, 2)
    if rng1 != rng
        error("ranges on lines 1 and 2 do not match.")
    end
    # parse line 3. Example: Variables: (a = 1, b = 2, c = 3)
    line = readline(io)
    m = match(r"Variables: (\([\w\s=,]+\))", line)
    if m === nothing
        error("expected Variables: at the start of line 3, got ", line, ".")
    end
    nt = let
        nt = parse_namedtuple(m.captures[1], 3)
        nms = collect(keys(nt))
        inds = collect(nt)
        if inds != 1:length(inds)
            if Set(inds) != Set(1:length(inds))
                error("indexes of variables on line 3 are not valid.")
            end
            tmp = Dict(i => n for (i, n) in zip(inds, nms))
            nt = NamedTuple{((tmp[i] for i = 1:length(inds))...,)}(1:length(inds))
        end
        nt
    end
    # parse line 4.  Example "(X) = Exogenous, (-) = Endogenous
    line = readline(io)
    m = match(r"\((.+?)\) = Exogenous, \((.+?)\) = Endogenous:", line)
    if m === nothing
        error("unexpected line 4: ", line, ".")
    end
    exog_mark, endo_mark = m.captures
    # parse line 5. Example "  NAME,  2000Q1:2010Q1, 2010Q2, 2010Q2:2020Q4"
    line = readline(io) * " "
    m = match(r"\s*NAME(.\S*)\s+([\w:]+(.\S*)\s*.*)", line)
    if m === nothing
        error("unexpected line 5: ", line, ".")
    end
    _name_delim = m.captures[1]
    _range_delim = m.captures[3]
    ranges = [parse_range(strip(str), true, 5) for str in split(m.captures[2], _range_delim; keepempty = false)]
    # parse the rest of it
    p = Plan{MIT{freq}}(rng, nt, falses(length(rng), length(nt)))
    ranges_inds = [[_offset(p, r)...] for r in ranges]
    pat = Regex("\\s+(\\w+)\\s*$(_name_delim)\\s*(.*)")
    for i = 1:length(nt)
        line = readline(io)
        m = match(pat, line)
        if m === nothing
            error("failed to parse line ", 5 + i, ": ", line)
        end
        var = Symbol(m.captures[1])
        var_ind = p.varshks[var]
        vals = split(m.captures[2], _range_delim; keepempty = false)
        for (rng_ind, val) in zip(ranges_inds, vals)
            val = strip(val)
            if val == exog_mark
                p.exogenous[rng_ind, var_ind] .= true
            elseif val == endo_mark
                nothing
            else
                error("unexpected $val on line ", 5 + i)
            end
        end
    end
    # @info " " freq rng nt ranges
    return p
end

function parse_range(str, allow_mit = true, line = nothing)
    e = Meta.parse(str)
    if occursin(":", str) || !allow_mit
        ans = Meta.isexpr(e, :call) && e.args[1] == :(:) ? eval(e) : nothing
        if !(ans isa UnitRange{<:MIT})
            error("expected range, got ", str, line === nothing ? "." : " on line $line.")
        end
    else
        ans = Meta.isexpr(e, :call) && e.args[1] == :(*) ? eval(e) : nothing
        if !(ans isa MIT)
            error("expected MIT, got ", str, line === nothing ? "." : " on line $line.")
        end
    end
    return ans
end

function parse_frequency(str, line = nothing)
    e = Meta.parse(str)
    ans = e isa Symbol ? eval(e) : nothing
    if !(ans isa Type && ans <: Frequency)
        error("expected frequency, got ", str, line === nothing ? "." : " on line $line.")
    end
    ans
end

function parse_namedtuple(str, line = nothing)
    e = Meta.parse(str)
    ans = Meta.isexpr(e, :tuple) && all(Meta.isexpr(ee, :(=)) for ee in e.args) ? eval(e) : nothing
    if !(ans isa NamedTuple{NAMES,NTuple{N,Int}} where {NAMES,N})
        error("expected NamedTuple, got ", str, line === nothing ? "." : " on line $line.")
    end
    ans
end

#######################################
# The user interface to modify plan instances.

"""
    setplanvalue!(plan, value, vars, range)
    setplanvalue!(plan, value, vars, date)
    setplanvalue!(plan, value, vars, dates)

Modify the status of the given variable(s) on the given date(s). If `value` is
`true` then variables become exogenous, otherwise they become endogenous.

"""
function setplanvalue!(p::Plan{T}, val::Bool, vars::Array{Symbol,1}, date::AbstractUnitRange{T}) where {T<:MIT}
    firstindex(p) ≤ first(date) && last(date) ≤ lastindex(p) || throw(BoundsError(p, date))
    idx1 = _offset(p, date)
    for v in vars
        idx2 = get(p.varshks, v, 0)
        if idx2 == 0
            throw(ArgumentError("Unknown variable or shock name: $v"))
        end
        p.exogenous[idx1, idx2] .= val
    end
    return p
end
setplanvalue!(p::Plan{MIT{Unit}}, val::Bool, vars::AbstractArray{Symbol,1}, date::AbstractUnitRange{Int}) = setplanvalue!(p, val, vars, UnitRange{MIT{Unit}}(date))
setplanvalue!(p::Plan{MIT{Unit}}, val::Bool, vars::AbstractArray{Symbol,1}, date::Integer) = setplanvalue!(p, val, vars, date*U:date*U)
setplanvalue!(p::Plan{MIT{Unit}}, val::Bool, vars::AbstractArray{Symbol,1}, date::MIT{Unit}) = setplanvalue!(p, val, vars, date:date)
setplanvalue!(p::Plan{T}, val::Bool, vars::AbstractArray{Symbol,1}, date::T) where {T<:MIT} = setplanvalue!(p, val, vars, date:date)
setplanvalue!(p::Plan, val::Bool, vars::AbstractArray{Symbol,1}, date) = (foreach(d -> setplanvalue!(p, val, vars, d), date); p)

"""
    exogenize!(plan, vars, date)

Modify the given plan so that the given variables will be exogenous on the given
dates. `vars` can be a `Symbol` or a `String` or a `Vector` of such. `date` can
be a moment in time (same type as the plan), or a range or an iterable or a
container.

"""
exogenize!(p::Plan, var, date) = setplanvalue!(p, true, Symbol[var,], date)
exogenize!(p::Plan, vars::AbstractVector, date) = setplanvalue!(p, true, Symbol[vars...], date)
exogenize!(p::Plan, vars::Tuple, date) = setplanvalue!(p, true, Symbol[vars...], date)

"""
    endogenize!(plan, vars, date)

Modify the given plan so that the given variables will be endogenous on the
given dates. `vars` can be a `Symbol` or a `String` or a `Vector` of such.
`date` can be a moment in time (same type as the plan), or a range or an
iterable or a container.

"""
endogenize!(p::Plan, var, date) = setplanvalue!(p, false, Symbol[var,], date)
endogenize!(p::Plan, vars::AbstractVector, date) = setplanvalue!(p, false, Symbol[vars...], date)
endogenize!(p::Plan, vars::Tuple, date) = setplanvalue!(p, false, Symbol[vars...], date)

"""
    exog_endo!(plan, exog_vars, endo_vars, date)

Modify the given plan so that the given variables listed in `exog_vars` will be
exogenous and the variables listed in `endo_vars` will be endogenous on the
given dates. `exog_vars` and `endo_vars` can each be a `Symbol` or a `String` or
a `Vector` of such. `date` can be a moment in time (same type as the plan), or a
range or an iterable or a container.

"""
function exog_endo!(p::Plan, exog, endo, date)
    setplanvalue!(p, true, Symbol[exog...], date)
    setplanvalue!(p, false, Symbol[endo...], date)
end
function exog_endo!(p::Plan, exog::Union{AbstractString,Symbol,ModelVariable}, endo::Union{AbstractString,Symbol,ModelVariable}, date)
    setplanvalue!(p, true, Symbol[exog], date)
    setplanvalue!(p, false, Symbol[endo], date)
end

"""
    endo_exog!(plan, endo_vars, exog_vars, date)

Modify the given plan so that the given variables listed in `exog_vars` will be
exogenous and the variables listed in `endo_vars` will be endogenous on the
given dates. `exog_vars` and `endo_vars` can each be a `Symbol` or a `String` or
a `Vector` of such. `date` can be a moment in time (same type as the plan), or a
range or an iterable or a container.

"""
@inline endo_exog!(p::Plan, endo, exog, date) = exog_endo!(p, exog, endo, date)

"""
autoexogenize!(plan, model, date)

Modify the given plan according to the "autoexogenize" protocol defined in the
given model. All variables in the autoexogenization list become endogenous and
their corresponding shocks become exogenous over the given date or range. `date`
can be a moment in time (same frequency as the given plan), a range, an
iterable, or a container.

"""
function autoexogenize!(p::Plan, m::Model, date)
    auto_vars = keys(m.autoexogenize)
    auto_shks = values(m.autoexogenize)
    exog_endo!(p, auto_vars, auto_shks, date)
end

#######################################
# The user interface to prepare data for simulation.

@inline TimeSeriesEcon.frequencyof(p::Plan) = frequencyof(p.range)
@inline TimeSeriesEcon.firstdate(p::Plan) = first(p.range)
@inline TimeSeriesEcon.lastdate(p::Plan) = last(p.range)

#######################################
# The internal interface to simulations code.

"""
    plansum(model, plan)

Return the total number of exogenous variables in the simulation plan. Periods
over which initial and final conditions are imposed are not counted in this sum.

"""
plansum(m::Model, p::Plan) = sum(p.exogenous[(1+m.maxlag):(end-m.maxlead), :])

export setexog!
"""
    setexog!(plan, t, vinds)

Modify the plan at time t such that `vinds` are exogenous and the rest are
endogenous.

"""
function setexog!(p::Plan, tt::Int, vinds)
    p.exogenous[tt, :] .= false
    p.exogenous[tt, vinds] .= true
end
setexog!(p::Plan{T}, tt::T, vinds) where {T<:MIT} = setexog!(p, _offset(p, tt), vinds)

end # module Plans

using .Plans
export Plan,
    exogenize!, endogenize!,
    exog_endo!, endo_exog!,
    autoexogenize!,
    exportplan, importplan


text:
module ExprTools
using MLStyle.MatchCore
export take_type_parameters!, get_type_parameters, get_type_parameters_ordered
export @reexport

@nospecialize
function take_type_parameters!(syms, ex)::Nothing
    @sswitch ex begin
        @case :($a >: $_) || :($a <: $_)
        @assert a isa Symbol
        push!(syms, a)
        return
        @case :($_ >: $b >: $_) || :($_ <: $b <: $_)
        @assert b isa Symbol
        push!(syms, b)
        return
        @case ::Symbol
        push!(syms, ex)
        return
        @case _
        return
    end
end

function get_type_parameters(args::AbstractArray{T, 1})::AbstractSet{Symbol} where {T}
    syms = Set{Symbol}()
    for arg in args
        take_type_parameters!(syms, arg)
    end
    syms
end

function get_type_parameters_ordered(args::AbstractArray{T, 1})::Vector{Symbol} where {T}
    syms = Symbol[]
    for arg in args
        take_type_parameters!(syms, arg)
    end
    unique!(syms)
    syms
end

macro reexport(m)
    m = __module__.eval(m)
    ns = names(m)
    m_name = nameof(m)
    ns = [n for n in ns if n !== m_name]
    isempty(ns) ? nothing : esc(:(export $(ns...)))
end
@specialize

end

ext: jl


text:
"""
    pedlist,idtable = read_ped(filename [, integer, order, hasupg])
    pedlist,idtable = read_ped(T, filename [, integer, order, hasupg])

Reads a pedigree file and store the pedigree list to an integer matrix `pedlist`
(2 x *n*; where *n* is the number of animals found in the pedigree). The argument
`T` defines the data type of `pedlist` and the default is Int. This function
supports two types of ID: character and integer.

By default, this function reads the original as characters and converts it to an
integer. The table is in a dictionary `idtable`, and it comes as the 2nd return-value.
If the ID is an integer, you can put the option `integer=true`, where the function
simply reads it as an integer without any conversion. No table is created so you will
have a single return value, `pedlist`.

With integer coding, the function treats all sires and dams not in the animal field
as real animals and creates new entries for them. It is not appropriate when your
pedigree file is "renumbered" by RENUMF90 with unknown parent groups (UPGs) where
UPG code is larger than the ID for real animals. The option `hasupg=true` does
nothing for UPG in sire/dam fields so that such UPG-code stays in `pedlist` as is.
The pedigree file has a space-separated text format with at least three fields: animal,
sire, and dam ID. By default, the first three fields are assigned to these IDs.

The optional argument `order` changes the position of fields. The default is
`order=[1,2,3]` and the 1st field for animal, the 2nd for sire, and the 3rd for dam.

### Examples
```jldoctest
julia> using PedigreeTools;

# for character data
julia> pedlist,idtable = read_ped("pedigree.txt");

# for character data with Int32 storage
julia> pedlist,idtable = read_ped(Int32,"pedigree.txt");

# for integer coding
julia> pedlist = read_ped("pedigree.txt", integer=true);
```
"""
function read_ped(Ti::DataType, filename::String; integer=false, order::Vector{To}=[1,2,3], hasupg=false) where {To<:Integer}
   if count(order.<1)>0
      error("`order` has 0 or negative value(s)")
   end
   if !(Ti<:Integer)
      error("non integer datatype specified")
   end
   if integer
      return read_ped_integer(Ti, filename, order=order; hasupg=hasupg)
   else
      return read_ped_general(Ti, filename, order=order)
   end
end

function read_ped(filename::String, Ti::DataType=Int; integer=false, order::Vector{To}=[1,2,3], hasupg=false) where {To<:Integer}
   return read_ped(Ti, filename, integer=integer, order=order, hasupg=hasupg)
end

# a general pedigree file; including characters
# change the name from ReadPedigreeGeneral to read_ped_general
function read_ped_general(Ti::DataType, filename::String; order::Vector{To}=[1,2,3]) where To<:Integer
   idtable = Dict{String,Ti}()
   idanim = Dict{String,Ti}()
   #pedtuple = Array{Tuple{Ti,Ti,Ti},1}()
   pedtuple = Array{NTuple{3,Ti},1}()
   line = ""
   pedline = zeros(Ti,3)
   lastid = zero(Ti)

   # checks
   if !isfile(filename)
      error("file not found")
   end

   # read the file
   idtable["0"] = 0
   open(filename,"r") do ios
      while(!eof(ios))
         line = strip(readline(ios))
         tokens = split(line)
         # animal
         id = tokens[order[1]]
         if haskey(idanim,id)
            @warn "ID: '"*id*"' appears in this file again. The last definision will be used."
            idanim[id] = idanim[id] + 1
         else
            idanim[id] = 1
         end
         # read the tokens
         for i=1:3
            id = tokens[order[i]]
            # no existing key
            if !haskey(idtable,id)
               lastid = lastid + 1
               idtable[id] = lastid
            end
            # set integer code
            pedline[i] = idtable[id]
         end
         # make a tuple
         push!(pedtuple,tuple(pedline...))
      end
   end

   # convert tuple to list
   pedlist = convert_pedigree_tuple_to_matrix(pedtuple,lastid)

   return pedlist,idtable
end

# an integer pedigree file; including characters
# change the name from ReadPedigreeInteger to read_ped_integer
function read_ped_integer(Ti::DataType, filename::String; order::Vector{To}=[1,2,3], hasupg=false) where To<:Integer
   #pedtuple = Array{Tuple{Ti,Ti,Ti},1}()
   pedtuple = Array{NTuple{3,Ti},1}()
   idanim = Dict{Ti,Ti}()
   lastid = zero(Ti)
   i = zero(Ti)
   k = zero(Ti)
   pedline = zeros(Ti,3)

   # checks
   if !isfile(filename)
      error("file not found")
   end

   # switch how to the maximum id
   find_maxid1(list::Vector{Ti},lastid) where Ti<:Integer = max(lastid,list[1])
   find_maxid2(list::Vector{Ti},lastid) where Ti<:Integer = max(lastid,maximum(list))
   if hasupg
      # maximum id only from the "animal" column
      find_maxid = find_maxid1
   else
      # maximum id only from all the columns
      find_maxid = find_maxid2
   end

   # read the file
   lastid = 0
   open(filename,"r") do ios
      while(!eof(ios))
         line = strip(readline(ios))
         tokens = split(line)
         # animal
         id = parse(Int,tokens[order[1]])
         if haskey(idanim,id)
            @warn "ID: '"*string(id)*"' appears in this file again. The last definision will be used."
            idanim[id] = idanim[id] + 1
         else
            idanim[id] = 1
         end
         # read the tokens
         for j=1:3
            try
               pedline[j] = parse(Int,tokens[order[j]])
            catch
               error("conversion error (non-integer value found): " * tokens[order[j]])
            end
         end
         # last id
         lastid = find_maxid(pedline,lastid)
         # make a tuple
         push!(pedtuple,tuple(pedline...))
      end
   end

   # convert tuple to list
   pedlist = convert_pedigree_tuple_to_matrix(pedtuple,lastid)

   return pedlist
end

# pedigree tuple to pedigree list
function convert_pedigree_tuple_to_matrix(pedtuple::Array{NTuple{3,Ti},1},lastid::Ti) where {Ti<:Integer}
   pedlist = zeros(Ti,2,lastid)
   for p in pedtuple
      if p[1]<1 || p[1]>lastid
         error("Invalid ID found")
      end
      pedlist[1:2,collect(p[1])] = collect(p[2:3])
   end
   return pedlist
end

# check pedigree
# change the name from CheckPedigree to check_ped
"""
    bool = check_ped(pedlist [, parentsfirst, warning])

Checks if the pedigree list has no apparent error. The pedigree list `pedlist` is
an integer matrix generated with `ReadPedigree` from a file. This function checks
the following errors.

* Sires (or dams) appear as dams (or sires).
* The list has a loop i.e., an animal's ancestor is also progeny of the animal.
* Parents precede their progeny i.e., IDs are chronologically sorted (optional; checked with `parentsfirst=true`).

It returns `true` with no error but `false` with a warning message where the error is (omit with `warining=false`).

### Examples

```jldoctest
julia> using PedigreeTools;

julia> pedlist,idtable = read_ped("pedigree.txt");

julia> check_ped(pedlist)
true

# checks if parents precede their progeny
julia> check_ped(pedlist, parentsfirst=true)
true
```
"""
function check_ped(pedlist::Matrix{Ti}; parentsfirst::Bool=false, warning::Bool=true) where Ti<:Integer
   n::Ti = size(pedlist,2)
   mark = zeros(Bool,n)

   # check sire/dam inconsistency
   mark = zeros(Bool,n)
   for i=1:n
      k = pedlist[1,i]
      if k>=1 && k<=n
         mark[k] = true
      end
   end
   for i=1:n
      k = pedlist[2,i]
      if k>=1 && k<=n
         if mark[k]
            if warning
               @warn "sire/dam found in the alternative sex: " * string(k)
            end
            return false
         end
      end
   end

   # check pedigree loop
   for i=1:n
      loop::Bool = hasloop(pedlist,i,i,0,false)
      if loop
         if warning
            @warn "pedigree loops found: " * string(i)
         end
         return false
      end
   end

   # check if sorted chronologically
   if parentsfirst
      for i=1:n
         unsorted::Bool = csorted(pedlist,i,i,false)
         if unsorted
            if warning
               @warn "pedigree not chronologically sorted: " * string(i)
            end
            return false
         end
      end
   end

   return true
end

# recursive tracing back of pedigree
function hasloop(pedlist::Matrix{Ti},ini::Ti,i::Ti,gen::Int,loop::Bool) where Ti<:Integer
   # immediate exit if it has a loop
   if loop==true
      return true
   end
   # existing animal
   if i>=1 && i<=size(pedlist,2)
      if gen>0 && i==ini
         # the same animal found as ancestor = pedigree loop
         return true
      else
         # tracing back
         return hasloop(pedlist,ini,pedlist[1,i],gen+1,loop) |
                hasloop(pedlist,ini,pedlist[2,i],gen+1,loop)
      end
   else
      # unknown base animal; do nothing
      return loop
   end
end

# check if chronologically sorted
function csorted(pedlist::Matrix{Ti},ini::Ti,i::Ti,unsorted::Bool) where Ti<:Integer
   # immediate exit if unsorted
   if unsorted==true
      return true
   end
   # existing animal
   if i>=1 && i<=size(pedlist,2)
      if i > ini
         # ancestor with larger ID than the staring animal = unsorted
         return true
      else
         # tracing back
         return csorted(pedlist,ini,pedlist[1,i],unsorted) |
                csorted(pedlist,ini,pedlist[2,i],unsorted)
      end
   else
      # unknown base animal; do nothing
      return unsorted
   end
end

"""
    newpedlist,newidtable = set_upg(pedlist, idtable, upgtable)

Assigns unknown parent groups (UPG) to unknown parents following the rule
in `upgtable` for a pedigree list `pedlist` and the original ID table `idtable`
See the following example;
5 UPG and the animals with IDs "A", "B", and "C" are grouped into the 1st UPG and "D" and "E" are grouped into the 2nd UPG.

```
upgtable = Dict("A"=>1, "B"=>1, "C"=>1, "D"=>2, "E"=>2)
```

The "animals" defined as UPG in `upgtable` are eliminated from the pedigree list.
If there are parent animals, they are replaced with integer UPG-code larger than the number of real animals. The table `idtable` is also updated accordingly.
"""
function set_upg(pedlist::Matrix{Ti},idtable::Dict{String,Ti},upgtable::Dict{String,Tj}; sortkeys=false) where {Ti<:Integer,Tj}
   #newpedtuple = Array{Tuple{Int,Int,Int},1}()
   newpedtuple = Array{NTuple{3,Ti},1}()
   invidtable = Dict{Ti,String}()
   newidtable = Dict{String,Ti}()
   upgdef = Dict{Tj,Ti}()
   pedline = zeros(Ti,3)
   lastcode = zero(Ti)
   # sort
   if sortkeys
      sortfunc(x) = sort(collect(x))
   else
      sortfunc = x->x
   end
   # the number of UPGs
   nupg=0
   for upgid in sortfunc(keys(upgtable))
      if upgid!="0" && haskey(idtable,upgid)
         nupg = nupg + 1
      end
   end
   # the number of real animals
   n = size(pedlist,2)
   m = n - nupg
   if m<1
      error("No real animals")
   end
   lastupg = m
   # make a new table with upg
   invidtable[0] = "0"
   newidtable["0"] = 0
   for id in sortfunc(keys(idtable))
      if id!="0"
         code = idtable[id]
         invidtable[code] = id
         if haskey(upgtable,id)
            # this id is upg: assign a large number (>m)
            if !haskey(newidtable,id)
               # grup ID
               group = upgtable[id]
               if haskey(upgdef,group)
                  # exsiting definition
                  newidtable[id] = upgdef[group]
               else
                  # new definition
                  lastupg = lastupg + 1
                  newidtable[id] = lastupg
                  upgdef[group] = lastupg
               end
            end
         else
            # this id is real
            # no existing key
            if !haskey(newidtable,id)
               lastcode = lastcode + 1
               newidtable[id] = lastcode
            end
         end
      end
   end
   # make a new pedigree list
   for id in sortfunc(keys(idtable))
      if id!="0"
         # animal code with id
         code = idtable[id]
         # store only if it is a real animal
         if !haskey(upgtable,id)
            # sire id and dam id
            sid = invidtable[pedlist[1,code]]
            did = invidtable[pedlist[2,code]]
            # new code
            pedline[1] = newidtable[id]
            pedline[2] = newidtable[sid]
            pedline[3] = newidtable[did]
         end
         # make a tuple
         push!(newpedtuple,tuple(pedline...))
      end
   end

   # convert tuple to list
   newpedlist = convert_pedigree_tuple_to_matrix(newpedtuple,lastcode)

   return newpedlist,newidtable
end


text:
export build_tarballs, autobuild, print_buildjl, product_hashes_from_github_release
import GitHub: gh_get_json, DEFAULT_API
import SHA: sha256

"""
    build_tarballs(ARGS, src_name, src_version, sources, script, platforms,
                   products, dependencies)

This should be the top-level function called from a `build_tarballs.jl` file.
It takes in the information baked into a `build_tarballs.jl` file such as the
`sources` to download, the `products` to build, etc... and will automatically
download, build and package the tarballs, generating a `build.jl` file when
appropriate.  Note that `ARGS` should be the top-level Julia `ARGS` command-
line arguments object.  This function does some rudimentary parsing of the
`ARGS`, call it with `--help` in the `ARGS` to see what it can do.
"""
function build_tarballs(ARGS, src_name, src_version, sources, script,
                        platforms, products, dependencies)
    # See if someone has passed in `--help`, and if so, give them the
    # assistance they so clearly long for
    if "--help" in ARGS
        println(strip("""
        Usage: build_tarballs.jl [target1,target2,...] [--only-buildjl]
                                 [--verbose] [--help]

        Options:
            targets         By default `build_tarballs.jl` will build a tarball
                            for every target within the `platforms` variable.
                            To override this, pass in a list of comma-separated
                            target triplets for each target to be built.  Note
                            that this can be used to build for platforms that
                            are not listed in the 'default list' of platforms
                            in the build_tarballs.jl script.

            --part=n/m      For breaking up long builds into multiple jobs,
                            divides the `platforms` list into `m` roughly
                            equal parts and then only builds part `n`
                            (`1 ≤ n ≤ m`).  (Does not produce `build.jl`,
                            which must be generated by a later
                            `--only-buildjl` stage.)

            --verbose       This streams compiler output to stdout during the
                            build which can be very helpful for finding bugs.
                            Note that it is colorized if you pass the
                            --color=yes option to julia, see examples below.

            --debug         This causes a failed build to drop into an
                            interactive bash shell for debugging purposes.

            --only-buildjl  This disables building of any tarballs, and merely
                            reconstructs a `build.jl` file from a github
                            release.  This is mostly useful as a later stage in
                            a travis/github releases autodeployment setup.

            --help          Print out this message.

        Examples:
            julia --color=yes build_tarballs.jl --verbose
                This builds all tarballs, with colorized output.

            julia build_tarballs.jl x86_64-linux-gnu,i686-linux-gnu
                This builds two tarballs for the two platforms given, with a
                minimum of output messages.
        """))
        return nothing
    end

    function check_flag(flag)
        flag_present = flag in ARGS
        ARGS = filter!(x -> x != flag, ARGS)
        return flag_present
    end

    # This sets whether we should build verbosely or not
    verbose = check_flag("--verbose")

    # This sets whether we drop into a debug shell on failure or not
    debug = check_flag("--debug")

    # This flag skips actually building and instead attempts to reconstruct a
    # build.jl from a GitHub release page.  Use this to automatically deploy a
    # build.jl file even when sharding targets across multiple CI builds.
    only_buildjl = check_flag("--only-buildjl")

    # --part=n/m builds only part n out of m divisions
    # of the platforms list.
    i = Compat.findlast(x -> startswith(x, "--part="), ARGS)
    should_override_platforms = i !== nothing
    if should_override_platforms
        if i != Compat.findfirst(x -> startswith(x, "--part="), ARGS)
            error("multiple --part arguments are not allowed")
        end
        p = parse.(Int, split(ARGS[i][8:end], '/'))
        (length(p) == 2 && p[2] > 0 && 1 ≤ p[1] ≤ p[2]) || error("invalid argument ", ARGS[i])
        n = (length(platforms) + p[2]-1) ÷ p[2]
        platforms = platforms[n*(p[1]-1)+1:min(end,n*p[1])]
        deleteat!(ARGS, i)
    end

    # If the user passed in a platform (or a few, comma-separated) on the
    # command-line, use that instead of our default platforms
    if length(ARGS) > 0
        should_override_platforms = true
        platforms = platform_key.(split(ARGS[1], ","))
    end

    # If we're running on CI (Travis, GitLab CI, etc...) and this is a
    # tagged release, automatically determine bin_path by building up a URL
    repo = get_repo_name()
    tag = get_tag_name()

    product_hashes = if !only_buildjl
        # If the user didn't just ask for a `build.jl`, go ahead and actually build
        Compat.@info("Building for $(join(triplet.(platforms), ", "))")

        # Build the given platforms using the given sources
        autobuild(pwd(), src_name, src_version, sources, script, platforms,
                         products, dependencies; verbose=verbose, debug=debug)
    else
        msg = strip("""
        Reconstructing product hashes from GitHub Release $(repo)/$(tag)
        """)
        Compat.@info(msg)

        # Reconstruct product_hashes from github
        product_hashes_from_github_release(repo, tag;
            product_filter="v$(src_version)",
            verbose=verbose
        )
    end

    # If we didn't override the default set of platforms OR we asked for only
    # a build.jl file, then write one out.  We don't write out when overriding
    # the default set of platforms because that is typically done either while
    # testing, or when we have sharded our tarball construction over multiple
    # invocations.
    if !should_override_platforms || only_buildjl
        # The location the binaries will be available from
        bin_path = "https://github.com/$(repo)/releases/download/$(tag)"

        # A dummy prefix to pass through products()
        dummy_products = products(Prefix(pwd()))
        print_buildjl(pwd(), src_name, src_version, dummy_products,
                      product_hashes, bin_path)

        if verbose
            Compat.@info("Writing out the following reconstructed build.jl:")
            print_buildjl(STDOUT, dummy_products, product_hashes, bin_path)
        end
    end

    return product_hashes
end

function build_tarballs(ARGS, src_name, sources, script, platforms, products,
                        dependencies)
    Compat.@warn("build_tarballs now requires a src_version parameter; assuming v\"1.0.0\"")
    return build_tarballs(
        ARGS,
        src_name,
        v"1.0.0",
        sources,
        script,
        platforms,
        products,
        dependencies
    )
end

# Helper function to get things from ENV, returning `nothing`
# if they either don't exist or are empty
function get_ENV(key)
    if !haskey(ENV, key)
        return nothing
    end

    if isempty(ENV[key])
        return nothing
    end

    return ENV[key]
end

function get_repo_name()
    # Helper function to synthesize repository slug from environment variables
    function get_gitlab_repo_name()
        owner = get_ENV("CI_REPO_OWNER")
        name = get_ENV("CI_REPO_NAME")
        if owner != nothing && name != nothing
            return "$(owner)/$(name)"
        end
        return nothing
    end

    # Helper function to guess repository slug from git remote URL
    function read_git_origin()
        try
            repo = LibGit2.GitRepo(".")
            url = LibGit2.url(LibGit2.get(LibGit2.GitRemote, repo, "origin"))
            owner = basename(dirname(url))
            if contains(owner, ":")
                owner = owner[findlast(owner, ':')+1:end]
            end
            name = basename(url)
            if endswith(name, ".git")
                name = name[1:end-4]
            end
            return "$(owner)/$(name)"
        catch
        end
        return nothing
    end

    return something(
        get_ENV("TRAVIS_REPO_SLUG"),
        get_gitlab_repo_name(),
        read_git_origin(),
        "<repo owner>/<repo name>",
    )
end

function get_tag_name()
    # Helper function to guess tag from current commit taggedness
    function read_git_tag()
        try
            repo = LibGit2.GitRepo(".")
            head_gitsha = LibGit2.GitHash(LibGit2.head(repo))
            for tag in LibGit2.tag_list(repo)
                tag_gitsha = LibGit2.GitHash(LibGit2.GitCommit(repo, tag))
                if head_gitsha == tag_gitsha
                    return tag
                end
            end
        catch
        end
        return nothing
    end

    return something(
        get_ENV("TRAVIS_TAG"),
        get_ENV("CI_COMMIT_TAG"),
        read_git_tag(),
        "<tag>",
    )
end

"""
    autobuild(dir::AbstractString, src_name::AbstractString,
              src_version::VersionNumber, sources::Vector,
              script::AbstractString, platforms::Vector,
              products::Function, dependencies::Vector;
              verbose::Bool = true, debug::Bool = false)

Runs the boiler plate code to download, build, and package a source package
for a list of platforms.  `src_name` represents the name of the source package
being built (and will set the name of the built tarballs), `platforms` is a
list of platforms to build for, `sources` is a list of tuples giving
`(url, hash)` of all sources to download and unpack before building begins,
`script` is a string representing a `bash` script to run to build the desired
products, which are listed as `Product` objects within the vector returned by
the `products` function. `dependencies` gives a list of dependencies that
provide `build.jl` files that should be installed before building begins to
allow this build process to depend on the results of another build process.
Setting `debug` to `true` will cause a failed build to drop into an interactive
shell so that the build can be inspected easily.
"""
function autobuild(dir::AbstractString,
                   src_name::AbstractString,
                   src_version::VersionNumber,
                   sources::Vector,
                   script::AbstractString,
                   platforms::Vector,
                   products::Function,
                   dependencies::Vector;
                   verbose::Bool = true, debug::Bool = false)
    # If we're on CI and we're not verbose, schedule a task to output a "." every few seconds
    if (haskey(ENV, "TRAVIS") || haskey(ENV, "CI")) && !verbose
        run_travis_busytask = true
        travis_busytask = @async begin
            # Don't let Travis think we're asleep...
            Compat.@info("Brewing a pot of coffee for Travis...")
            while run_travis_busytask
                sleep(4)
                print(".")
            end
        end
    end

    # This is what we'll eventually return
    product_hashes = Dict()

    # If we end up packaging any local directories into tarballs, we'll store them here
    mktempdir() do tempdir
        # First, download the source(s), store in ./downloads/
        downloads_dir = joinpath(dir, "downloads")
        try mkpath(downloads_dir) catch; end

        # We must prepare our sources.  Download them, hash them, etc...
        sources = Any[s for s in sources]
        for idx in 1:length(sources)
            # If the given source is a local path that is a directory, package it up and insert it into our sources
            if typeof(sources[idx]) <: AbstractString
                if !isdir(sources[idx])
                    error("Sources must either be a pair (url => hash) or a local directory")
                end

                # Package up this directory and calculate its hash
                tarball_path = joinpath(tempdir, basename(sources[idx]) * ".tar.gz")
                package(sources[idx], tarball_path; verbose=verbose)
                tarball_hash = open(tarball_path, "r") do f
                    bytes2hex(sha256(f))
                end

                # Now that it's packaged, store this into sources[idx]
                sources[idx] = (tarball_path => tarball_hash)
            elseif typeof(sources[idx]) <: Pair
                src_url, src_hash = sources[idx]

                # If it's a .git url, clone it
                if endswith(src_url, ".git")
                    src_path = joinpath(downloads_dir, basename(src_url))
                    if !isdir(src_path)
                        repo = LibGit2.clone(src_url, src_path; isbare=true)
                    else
                        LibGit2.with(LibGit2.GitRepo(src_path)) do repo
                            LibGit2.fetch(repo)
                        end
                    end
                else
                    if isfile(src_url)
                        # Immediately abspath() a src_url so we don't lose track of
                        # sources given to us with a relative path
                        src_path = abspath(src_url)

                        # And if this is a locally-sourced tarball, just verify
                        verify(src_path, src_hash; verbose=verbose)
                    else
                        # Otherwise, download and verify
                        src_path = joinpath(downloads_dir, basename(src_url))
                        download_verify(src_url, src_hash, src_path; verbose=verbose)
                    end
                end

                # Now that it's downloaded, store this into sources[idx]
                sources[idx] = (src_path => src_hash)
            else
                error("Sources must be either a `URL => hash` pair, or a path to a local directory")
            end
        end

        # Our build products will go into ./products
        out_path = joinpath(dir, "products")
        try mkpath(out_path) catch; end

        for platform in platforms
            target = triplet(platform)

            # We build in a platform-specific directory
            build_path = joinpath(pwd(), "build", target)
            try mkpath(build_path) catch; end

            cd(build_path) do
                src_paths, src_hashes = collect(zip(sources...))

                # Convert from tuples to arrays, if need be
                src_paths = collect(src_paths)
                src_hashes = collect(src_hashes)
                prefix, ur = setup_workspace(
                    build_path,
                    src_paths,
                    src_hashes,
                    dependencies,
                    platform;
                    verbose=verbose,
                    downloads_dir=downloads_dir
                )

                # Don't keep the downloads directory around
                rm(joinpath(prefix, "downloads"); force=true, recursive=true)

                # Collect dependency manifests so that our auditing doesn't touch these files that
                # were installed by dependencies
                manifest_dir = joinpath(prefix, "manifests")
                dep_manifests = if isdir(manifest_dir)
                   [joinpath(prefix, "manifests", f) for f in readdir(manifest_dir)]
                else
                    String[]
                end

                dep = Dependency(src_name, products(prefix), script, platform, prefix)
                build(ur, dep; verbose=verbose, autofix=true, ignore_manifests=dep_manifests, debug=debug)

                # Remove the files of any dependencies
                for dependency in dependencies
                    dep_script = script_for_dep(dependency)
                    m = Module(:__anon__)
                    eval(m, quote
                        using BinaryProvider
                        # Override BinaryProvider functionality so that it doesn't actually install anything
                        platform_key() = $platform
                        function write_deps_file(args...); end
                        function install(args...; kwargs...); end

                        # Include build.jl file to extract download_info
                        ARGS = [$(prefix.path)]
                        include_string($(dep_script))

                        # Grab the information we need in order to extract a manifest, then uninstall it
                        url, hash = download_info[platform_key()]
                        manifest_path = BinaryProvider.manifest_from_url(url; prefix=prefix)
                        BinaryProvider.uninstall(manifest_path; verbose=$verbose)
                    end)
                end

                # Once we're built up, go ahead and package this prefix out
                tarball_path, tarball_hash = package(
                    prefix,
                    joinpath(out_path, src_name),
                    src_version;
                    platform=platform,
                    verbose=verbose,
                    force=true,
               )
                product_hashes[target] = (basename(tarball_path), tarball_hash)

                # Destroy the workspace
                rm(dirname(prefix.path); recursive=true)
            end

            # If the whole build_path is empty, then remove it too.  If it's not, it's probably
            # because some other build is doing something simultaneously with this target, and we
            # don't want to mess with their stuff.
            if isempty(readdir(build_path))
                rm(build_path; recursive=true)
            end

            # Clean up this particular shard, so as not to run out of loopback devices
            unmount_shard(shards_dir(triplet(platform)))
        end
    end

    # At the end of all things, unmount all our shards so as to play nice with others
    unmount_all_shards()

    if (haskey(ENV, "TRAVIS") || haskey(ENV, "CI")) && !verbose
        run_travis_busytask = false
        wait(travis_busytask)
        println()
    end

    # Return our product hashes
    return product_hashes
end

function print_buildjl(io::IO, products::Vector, product_hashes::Dict,
                       bin_path::AbstractString)
    print(io, """
    using BinaryProvider # requires BinaryProvider 0.3.0 or later

    # Parse some basic command-line arguments
    const verbose = "--verbose" in ARGS
    const prefix = Prefix(get([a for a in ARGS if a != "--verbose"], 1, joinpath(@__DIR__, "usr")))
    """)

    # Print out products
    print(io, "products = [\n")
    for prod in products
        print(io, "    $(repr(prod)),\n")
    end
    print(io, "]\n\n")

    # Print binary locations/tarball hashes
    print(io, """
    # Download binaries from hosted location
    bin_prefix = "$bin_path"

    # Listing of files generated by BinaryBuilder:
    """)

    println(io, "download_info = Dict(")
    for platform in sort(collect(keys(product_hashes)))
        fname, hash = product_hashes[platform]
        pkey = platform_key(platform)
        println(io, "    $(pkey) => (\"\$bin_prefix/$(fname)\", \"$(hash)\"),")
    end
    println(io, ")\n")

    print(io, """
    # Install unsatisfied or updated dependencies:
    unsatisfied = any(!satisfied(p; verbose=verbose) for p in products)
    if haskey(download_info, platform_key())
        url, tarball_hash = download_info[platform_key()]
        if unsatisfied || !isinstalled(url, tarball_hash; prefix=prefix)
            # Download and install binaries
            install(url, tarball_hash; prefix=prefix, force=true, verbose=verbose)
        end
    elseif unsatisfied
        # If we don't have a BinaryProvider-compatible .tar.gz to download, complain.
        # Alternatively, you could attempt to install from a separate provider,
        # build from source or something even more ambitious here.
        error("Your platform \$(triplet(platform_key())) is not supported by this package!")
    end

    # Write out a deps.jl file that will contain mappings for our products
    write_deps_file(joinpath(@__DIR__, "deps.jl"), products)
    """)
end

function print_buildjl(build_dir::AbstractString, src_name::AbstractString,
                       src_version::VersionNumber, products::Vector,
                       product_hashes::Dict, bin_path::AbstractString)
    mkpath(joinpath(build_dir, "products"))
    open(joinpath(build_dir, "products", "build_$(src_name).v$(src_version).jl"), "w") do io
        print_buildjl(io, products, product_hashes, bin_path)
    end
end

"""
If you have a sharded build on Github, it would be nice if we could get an auto-generated
`build.jl` just like if we build serially.  This function eases the pain by reconstructing
it from a releases page.
"""
function product_hashes_from_github_release(repo_name::AbstractString, tag_name::AbstractString;
                                            product_filter::AbstractString = "",
                                            verbose::Bool = false)
    # Get list of files within this release
    release = gh_get_json(DEFAULT_API, "/repos/$(repo_name)/releases/tags/$(tag_name)", auth=github_auth())

    # Try to extract the platform key from each, use that to find all tarballs
    function can_extract_platform(filename)
        # Short-circuit build.jl because that's quite often there.  :P
        if startswith(filename, "build") && endswith(filename, ".jl")
            return false
        end

        unknown_platform = typeof(extract_platform_key(filename)) <: UnknownPlatform
        if unknown_platform && verbose
            Compat.@info("Ignoring file $(filename); can't extract its platform key")
        end
        return !unknown_platform
    end
    assets = [a for a in release["assets"] if can_extract_platform(a["name"])]
    assets = [a for a in assets if contains(a["name"], product_filter)]

    # Download each tarball, hash it, and reconstruct product_hashes.
    product_hashes = Dict()
    mktempdir() do d
        for asset in assets
            # For each asset (tarball), download it
            filepath = joinpath(d, asset["name"])
            url = asset["browser_download_url"]
            BinaryProvider.download(url, filepath; verbose=verbose)

            # Hash it
            hash = open(filepath) do file
                return bytes2hex(sha256(file))
            end

            # Then fit it into our product_hashes
            file_triplet = triplet(extract_platform_key(asset["name"]))
            product_hashes[file_triplet] = (asset["name"], hash)

            if verbose
                Compat.@info("Calculated $hash for $(asset["name"])")
            end
        end
    end

    return product_hashes
end


text:
# Helpers for working with GloVe in Julia
module GloVe

using SparseArrays
using LinearAlgebra

export CREC, WORD_INFO

struct CREC
    word1::Int32
    word2::Int32
    val::Float64
end

struct WORD_INFO
    index::Int64
    count::Int64
end

# Load a vocabulary file as a Dict
function load_vocab(filename)
    str2idx = Dict{String, WORD_INFO}()
    idx2str = Array{String, 1}()
    open(filename) do f
        for (i, line) in enumerate(eachline(f))
            entry = split(line)
            str2idx[entry[1]] = WORD_INFO(i, parse(Int64, entry[2]))
            push!(idx2str, entry[1])
        end
    end
    return (str2idx, idx2str)
end

# Load the GloVe binary vectors
function load_bin_vectors(filename, V)
    n = div(filesize(filename), 8) # Total number of params (8 bytes per real)
    dim = div(n - 2*V, 2*V) # Solve for dimension of embedding
    W = zeros(V, dim)
    U = zeros(V, dim)
    b_w = zeros(V)
    b_u = zeros(V)
    open(filename) do f
        # Read Word Vectors and Biases
        for i = 1:V
            for j = 1:dim
                W[i, j] = read(f, Float64)
            end
            b_w[i] = read(f, Float64)
        end
        # Read context Vectors and Biases
        for i = 1:V
            for j = 1:dim
                U[i, j] = read(f, Float64)
            end
            b_u[i] = read(f, Float64)
        end
    end
    return (W, b_w, U, b_u)
end

# Convenience wrapper to simplify loading
function load_model(embedding_path)
    vocab_path = "$(dirname(embedding_path))/vocab$(match(r"-C[0-9]+-V[0-9]+", embedding_path).match).txt"
    vocab, ivocab = load_vocab(vocab_path)
    d = parse(Int64, match(r"-W[0-9]+", embedding_path).match[3:end])  # window
    V = length(vocab)
    W, b_w, U, b_u = load_bin_vectors(embedding_path, V)
    D = size(W, 2)
    return (vocab=vocab, ivocab=ivocab, W=W, b_w=b_w, U=U, b_u=b_u, V=V, D=D, d=d,
        vocab_path=abspath(vocab_path), embedding_path=abspath(embedding_path))
end

# Convenience wrapper to simplify loading
function load_model(embedding_path, vocab_path)
    vocab, ivocab = load_vocab(vocab_path)
    d = parse(Int64, match(r"-W[0-9]+", embedding_path).match[3:end])  # window
    V = length(vocab)
    W, b_w, U, b_u = load_bin_vectors(embedding_path, V)
    D = size(W, 2)
    return (vocab=vocab, ivocab=ivocab, W=W, b_w=b_w, U=U, b_u=b_u, V=V, D=D, d=d,
        vocab_path=abspath(vocab_path), embedding_path=abspath(embedding_path))
end

# Save vectors in a text format (useful for running annalogy tests)
function save_text_vectors(filename, W, idx2str)
    (V, D) = size(W)
    open(filename, "w") do f
        for i = 1:V
            print(f, idx2str[i])
            for j = 1:D
                print(f, " ")
                print(f, W[i, j])
            end
            print(f, "\n")
        end
    end
end

# Save vectors in a binary format
function save_bin_vectors(filename, W, b_w, U, b_u)
    (V, dim) = size(W)
    open(filename, "w") do f
        # Read Word Vectors and Biases
        for i = 1:V
            for j = 1:dim
                write(f, W[i, j])
            end
            write(f, b_w[i])
        end
        # Read context Vectors and Biases
        for i = 1:V
            for j = 1:dim
                write(f, U[i, j])
            end
            write(f, b_u[i])
        end
    end
end

# Read a coocurence from file
function read_cooc(io)::Tuple{Int64,Int64,Float64}
    i::Int64 = read(io, Int32)
    j::Int64 = read(io, Int32)
    x::Float64 = read(io, Float64)
    return (i, j, x)
end

# Load the full coc_matrix.
function load_cooc(filename, V)::SparseMatrixCSC{Float64, Int64}
    I = Array{Int64,1}()
    J = Array{Int64,1}()
    X = Array{Float64,1}()
    open(filename) do f
        while (!eof(f))
            (i, j, x) = read_cooc(f)
            push!(I, i); push!(J, j); push!(X, x)
        end
    end
    return sparse(I, J, X, V, V)
end

# Load only specific rows/columns of the coocurence matrix
function load_cooc(filename, V, target_inds::Array{Int64,1})::SparseMatrixCSC{Float64, Int64}
    I = Array{Int64,1}()
    J = Array{Int64,1}()
    X = Array{Float64,1}()
    open(filename) do f
        while (!eof(f))
            (i, j, x) = read_cooc(f)
            if (i in target_inds || j in target_inds)
                push!(I, i); push!(J, j); push!(X, x)
            end
        end
    end
    return sparse(I, J, X, V, V)
end

# Replicated the GloVe preprocessing
function parse_coocs(text::String, vocab::Dict{String,GloVe.WORD_INFO}, window::Int64)::Tuple{Array{Int64,1},Array{Int64,1},Array{Float64,1}}
    words = split(text)
    i::Int64 = -1  # row index
    j::Int64 = -1  # col index
    l1::Int64 = 1  # position of center word
    l2::Int64 = 0  # position of context word
    net_offset::Int64 = 0  # offset (l1 - l2), excluding out-of-vocab "gaps"
    I = Array{Int64,1}()
    J = Array{Int64,1}()
    vals = Array{Float64,1}()
    # slide the center position
    for (l1, word) in enumerate(words)
        i = get(vocab, word, WORD_INFO(-1, -1)).index
        i == -1 && continue  # skip position if out-of-vocab
        l2 = l1  # align
        net_offset = 0  # reset
        while net_offset < window
            l2 -= 1  # increment
            l2 <= 0 && break
            j = get(vocab, words[l2], WORD_INFO(-1, -1)).index
            j == -1 && continue
            net_offset += 1  # word in-vocab, so increment net offset
            push!(I, i);
            push!(J, j);
            push!(vals, 1.0/net_offset)
        end
    end
    return vcat(I, J), vcat(J, I), vcat(vals, vals)
end

# Build a sparse cooc matrix from a single document
function doc2cooc(text::String, vocab::Dict{String,GloVe.WORD_INFO}, window::Int64)::SparseMatrixCSC{Float64,Int64}
    V = length(vocab)
    I, J, vals = parse_coocs(text, vocab, window)
    return sparse(I, J, vals, V, V)
end

# Build a sparse cooc matrix from an iterable set of documents
function docs2cooc(texts, vocab::Dict{String,GloVe.WORD_INFO}, window::Int64;verbose=false)::SparseMatrixCSC{Float64,Int64}
    V = length(vocab)
    I = Array{Int64,1}()
    J = Array{Int64,1}()
    vals = Array{Float64,1}()
    for (num, text) in enumerate(texts)
        i, j, val = parse_coocs(text, vocab, window)
        append!(I, i); append!(J, j); append!(vals, val)
        if (verbose && num % 1000 == 0)
            println("Parsing document $num")
        end
    end
    return sparse(I, J, vals, V, V)
end

# Save a coocurence matrix to file
function save_coocs(out_file::String, X::SparseMatrixCSC{Float64, Int64})
    open(out_file, "w") do out_io
        for (i, j, x) in zip(findnz(X)...)
            write(out_io, Int32(i))
            write(out_io, Int32(j))
            write(out_io, x)
        end
    end
end

# Stream co-occurrences, adding the perturbations
function perturb_coocs(in_file::String, out_file::String,
        X::SparseMatrixCSC{Float64, Int64})::Int64
    local n_perts = 0
    local i::Int32
    local j::Int32
    local x::Float64
    local Xij::Float64
    open(in_file) do in_io
        open(out_file, "w") do out_io
            while !eof(in_io)
                i = read(in_io, Int32); j = read(in_io, Int32); x = read(in_io, Float64)
                Xij = X[i,j]
                if (Xij != 0.0)
                    n_perts += 1
                    if (x + Xij > 0.0) # Don't write zeros/neg; causes error
                        write(out_io, i)
                        write(out_io, j)
                        write(out_io, x + Xij)
                    end
                else
                    write(out_io, i)
                    write(out_io, j)
                    write(out_io, x)
                end
            end
        end
    end
    return n_perts
end

# GloVe weighting function
function f(x::Float64; max::Float64=100.0, alpha::Float64=0.75)::Float64
    if x >= max
        return 1.0
    else
        return (x/max)^alpha
    end
end

# Pointwise loss
function Li(W::AbstractArray, b_w::Array{Float64,1}, U::Array{Float64,2},
            b_u::Array{Float64,1}, X::SparseMatrixCSC{Float64}, i::Int64)
    (V, D) = size(U)
    Xi = X[i, :].nzval
    Ji = X[i, :].nzind
    if length(Xi) == 0
        return 0.0
    end
    diff = U[Ji, :] * W[i, :] + b_u[Ji] .+ b_w[i] - log.(Xi)
    return sum(f.(Xi) .* diff.^2)
end

# Total loss
function J(W::AbstractArray, b_w::Array{Float64,1}, U::Array{Float64,2},
            b_u::Array{Float64,1}, X::SparseMatrixCSC{Float64})
    (V, D) = size(U)
    tot = 0.0
    for i in 1:V
        tot = tot + Li(W, b_w, U, b_u, X, i)
    end
    return tot
end

# Gradient of the pointwise loss w.r.t. Wi
function ∇Li(W::Array{Float64,2}, b_w::Array{Float64,1},
        U::Array{Float64,2}, b_u::Array{Float64,1}, X::SparseMatrixCSC{Float64,Int64}, i::Int64)::Array{Float64,1}
    Xi = X[i, :].nzval; Ji = X[i, :].nzind;
    diff = U[Ji, :] * W[i, :] + b_u[Ji] .+ b_w[i] - log.(Xi)
    return ((2 * GloVe.f.(Xi) .* diff)' * U[Ji, :])'
end

# Hessian of the pointwise loss w.r.t. Wi
function ∇²Li(U::Array{Float64,2}, X::SparseMatrixCSC{Float64,Int64},
        i::Int64)::Array{Float64,2}
    Xi = X[i, :].nzval; Ji = X[i, :].nzind;
    temp = sqrt.(2.0 * GloVe.f.(Xi)) .* @view U[Ji, :]
    return temp' * temp
end

# Compute change in embedding due to document removal
function compute_IF_deltas(document::String, M, X::SparseMatrixCSC)
    Y = GloVe.doc2cooc(document, M.vocab, M.d)  # The perturbation
    affected_inds = unique(Y.rowval)  # Set of affected word vectors
    N = length(affected_inds)  # Number of vectors that will change
    deltas = Dict{Int64,Array{Float64,1}}()  # Hash deltas by word index
    if N != 0
        X̃ = dropzeros(max.(0.0, X - Y))  # non-neg catch incase doc out of corpus
        for idx in affected_inds
            gi = ∇Li(M.W, M.b_w, M.U, M.b_u, X, idx)  # original gradient
            Hi = ∇²Li(M.U, X, idx)  # original hessian
            g̃i = ∇Li(M.W, M.b_w, M.U, M.b_u, X̃, idx)  # perturbed gradient
            deltas[idx] = inv(Hi) * (gi - g̃i)
        end
    end
    return deltas
end

# Speed-up: avoid computing anything you don't need to
function compute_IF_deltas(document::String, M, X::SparseMatrixCSC,
        target_indices::Array{Int64,1}, inv_hessians::Dict{Int64,Array{Float64,2}},
        gradients::Dict{Int64,Array{Float64,1}})
    Y = GloVe.doc2cooc(document, M.vocab, M.d)
    select_affected_inds = intersect(unique(Y.rowval), target_indices)  # consider less
    N = length(select_affected_inds)  # Number of vectors that will change
    deltas = Dict{Int64,Array{Float64,1}}() # Deltas for each word vector
    if N != 0
        X̃ = dropzeros(max.(0.0, X - Y))  # non-neg catch incase doc out of corpus
        for idx in select_affected_inds
            g̃i = ∇Li(M.W, M.b_w, M.U, M.b_u, X̃, idx)
            deltas[idx] = inv_hessians[idx] * (gradients[idx] - g̃i)
        end
    end
    return deltas
end

# Where we've already loaded the perturbation as a cooc matrix
function compute_IF_deltas(Y::SparseMatrixCSC, M, X::SparseMatrixCSC,
        target_indices::Array{Int64,1}, inv_hessians::Dict{Int64,Array{Float64,2}},
        gradients::Dict{Int64,Array{Float64,1}})
    select_affected_inds = intersect(unique(Y.rowval), target_indices)  # consider less
    N = length(select_affected_inds)  # Number of vectors that will change
    deltas = Dict{Int64,Array{Float64,1}}() # Deltas for each word vector
    if N != 0
        X̃ = dropzeros(max.(0.0, X - Y))  # non-neg catch incase doc out of corpus
        for idx in select_affected_inds
            g̃i = ∇Li(M.W, M.b_w, M.U, M.b_u, X̃, idx)
            deltas[idx] = inv_hessians[idx] * (gradients[idx] - g̃i)
        end
    end
    return deltas
end

# Precompute a batch of gradients
function gradients_for(target_indices, M, X)
    num_words = length(target_indices)
    G = Dict{Int64,Array{Float64,1}}()
    for i in target_indices
        G[i] = ∇Li(M.W, M.b_w, M.U, M.b_u, X, i)
    end
    return G
end

# Precompute a batch of inverse hessians
function inv_hessians_for(target_indices, M, X)
    num_words = length(target_indices)
    H = Dict{Int64,Array{Float64,2}}()
    for i in target_indices
        H[i] = inv(∇²Li(M.U, X, i))
    end
    return H
end

end


text:
using .TestUtils

@testgroup "runtests" begin

constant_fixture1 = 1:2
constant_fixture2 = ["a", "b"]

global_fixture1_setup = []
global_fixture1_torndown = []
global_fixture1_vals = [10, 20]
global_fixture1 = @global_fixture for val in global_fixture1_vals
    @assert !(val in global_fixture1_setup)
    push!(global_fixture1_setup, val)
    @produce val
    push!(global_fixture1_torndown, val)
end

global_fixture2_setup = []
global_fixture2_torndown = []
global_fixture2_vals = ["x", "y"]
global_fixture2 = @global_fixture for val in global_fixture2_vals
    @assert !(val in global_fixture2_setup)
    push!(global_fixture2_setup, val)
    @produce val
    push!(global_fixture2_torndown, val)
end

# Constant + constant

c12_results = []
@testcase "c12" for x1 in constant_fixture1, x2 in constant_fixture2
    push!(c12_results, (x1, x2))
end

@testcase "check c12" begin
    @test c12_results == collect(Jute.rowmajor_product(constant_fixture1, constant_fixture2))
end

# Only using global fixture 1

g1_results = []

@testcase "check g1 setup" begin
    @test global_fixture1_setup == []
    @test global_fixture1_torndown == []
end

@testcase "g1" for x1 in global_fixture1
    @test global_fixture1_setup == global_fixture1_vals
    push!(g1_results, x1)
end

@testcase "check g1" begin
    @test g1_results == collect(global_fixture1_vals)
end

@testcase "check g1 is still there" begin
    @test global_fixture1_torndown == []
end

# Using global fixtures 1 and 2

g12_results = []

@testcase "check g2 setup" begin
    @test global_fixture2_setup == []
    @test global_fixture2_torndown == []
end

@testcase "g12" for x1 in global_fixture1, x2 in global_fixture2
    @test global_fixture1_setup == global_fixture1_vals
    @test global_fixture2_setup == global_fixture2_vals
    push!(g12_results, (x1, x2))
end

@testcase "check g1 destroyed" begin
    @test global_fixture1_torndown == global_fixture1_vals
end

@testcase "check g2 is still there" begin
    @test global_fixture2_torndown == []
end

@testcase "check g12" begin
    @test g12_results == collect(Jute.rowmajor_product(global_fixture1_vals, global_fixture2_vals))
end

# Using only global fixture 2
g2_results = []

@testcase "g2" for x2 in global_fixture2
    @test global_fixture2_setup == global_fixture2_vals
    push!(g2_results, x2)
end

@testcase "check g2" begin
    @test g2_results == collect(global_fixture2_vals)
end

@testcase "check g2 is destroyed" begin
    @test global_fixture2_torndown == global_fixture2_vals
end

# Check dependent global fixtures

gfs_state = Dict(key => 0 for key in ["a", "b", "c", "d"])

as = ["a1", "a2"]
bs = ["b1", "b2"]
cs = ["c1", "c2"]
combine_ab(a, b) = a * b
combine_ac(a, c) = a * c
combine_bc(b, c) = b * "+" * c

gf_as = @global_fixture for a in as
    gfs_state["a"] += 1
    @produce a
    gfs_state["a"] -= 1
end

gf_bs = @global_fixture for a in gf_as, b in bs
    total_values = 4
    @assert gfs_state["b"] >= 0 && gfs_state["b"] <= total_values - 1
    gfs_state["b"] += 1
    @produce combine_ab(a, b)
    @assert gfs_state["b"] >= 1 && gfs_state["b"] <= total_values
    gfs_state["b"] -= 1
end

gf_cs = @global_fixture for a in gf_as, c in cs
    total_values = 4
    @assert gfs_state["c"] >= 0 && gfs_state["c"] <= total_values - 1
    gfs_state["c"] += 1
    @produce combine_ac(a, c)
    @assert gfs_state["c"] >= 1 && gfs_state["c"] <= total_values
    gfs_state["c"] -= 1
end

gf_ds = @global_fixture for b in gf_bs, c in gf_cs
    total_values = 16
    @assert gfs_state["d"] >= 0 && gfs_state["d"] <= total_values - 1
    gfs_state["d"] += 1
    @produce combine_bc(b, c)
    @assert gfs_state["d"] >= 1 && gfs_state["d"] <= total_values
    gfs_state["d"] -= 1
end

combination_results = Array{Tuple{String, String}, 1}()

@testcase "combination" for d in gf_ds, b in gf_bs
    push!(combination_results, (d, b))
end

@testcase "check combination" begin
    cbs = [combine_ab(a, b) for (a, b) in Jute.rowmajor_product(as, bs)]
    ccs = [combine_ac(a, c) for (a, c) in Jute.rowmajor_product(as, cs)]
    cds = [combine_bc(b, c) for (b, c) in Jute.rowmajor_product(cbs, ccs)]
    @test combination_results == collect(Jute.rowmajor_product(cds, cbs))
end

@testcase "check all cleaned" begin
    @test all([gfs_state[key] == 0 for key in keys(gfs_state)])
end

# Check local fixtures

lf_sequence = []

lf_nodeps = @local_fixture begin
    push!(lf_sequence, "setup")
    @produce 1
    push!(lf_sequence, "teardown")
end

@testcase "lf nodeps test" for x in lf_nodeps
    push!(lf_sequence, "testcase $x")
end

@testcase "lf nodeps check" begin
    @test lf_sequence == ["setup", "testcase 1", "teardown"]
end

lf_sequence2 = []

gf_for_lf = @global_fixture for val in [1, 2]
    push!(lf_sequence2, "gf setup $val")
    @produce val
    push!(lf_sequence2, "gf teardown $val")
end

lf_nodeps = @local_fixture begin
    push!(lf_sequence2, "lf_nodeps setup")
    @produce "a"
    push!(lf_sequence2, "lf_nodeps teardown")
end

lf_deps = @local_fixture for x in lf_nodeps, y in 3:4, z in gf_for_lf
    push!(lf_sequence2, "lf_deps $x $y $z setup")
    @produce (x, y, z)
    push!(lf_sequence2, "lf_deps $x $y $z teardown")
end

@testcase "lf deps test" for x in lf_deps, y in lf_nodeps
    push!(lf_sequence2, "testcase $x $y")
end

@testcase "lf deps check" begin
    lf_sequence_ref = []
    push!(lf_sequence_ref, "gf setup 1")
    push!(lf_sequence_ref, "gf setup 2")
    for i in 3:4
        for gf in 1:2
            push!(lf_sequence_ref, "lf_nodeps setup")
            push!(lf_sequence_ref, "lf_deps a $i $gf setup")
            push!(lf_sequence_ref, "lf_nodeps setup")
            push!(lf_sequence_ref, "testcase $(("a", i, gf)) a")
            push!(lf_sequence_ref, "lf_nodeps teardown")
            push!(lf_sequence_ref, "lf_deps a $i $gf teardown")
            push!(lf_sequence_ref, "lf_nodeps teardown")
        end
    end
    push!(lf_sequence_ref, "gf teardown 1")
    push!(lf_sequence_ref, "gf teardown 2")

    @test lf_sequence2 == lf_sequence_ref
end

# Checks that the max_fails option works.
# Also checks that all remaining teardowns are called
# if runtests exists prematurely because max_fails was reached.
@testcase "test max fails" begin

    teardown_called = false
    tc1_executed = false
    tc3_executed = false

    gfx = @global_fixture begin
        @produce 1
        teardown_called = true
    end

    tcs = Jute.collect_testobjs() do
        @testcase "tc1" for x in gfx
            tc1_executed = true
        end

        @testcase "tc2" begin
            @test 1 == 2
        end

        @testcase "tc3" for x in gfx
            tc3_executed = true
        end
    end

    exitcode, output = nested_run_with_output(tcs, Dict(:max_fails => 1))
    @test tc1_executed
    @test teardown_called
    @test !tc3_executed
    @test exitcode == 1
end

# Checks that the instant teardown global fixture
# is actually torn down right after instantiation.
@testcase "test instant teardown" begin
    teardown_called = false
    tc1_executed = false

    gfx = @global_fixture instant_teardown=true begin
        @produce 1
        teardown_called = true
    end

    tcs = Jute.collect_testobjs() do
        @testcase "tc1" for x in gfx
            @test teardown_called
            tc1_executed = true
        end
    end

    exitcode = nested_run(tcs)
    @test teardown_called
    @test tc1_executed
    @test exitcode == 0
end

# Checks that the case of no testcases to run is handled appropriately.
@testcase "test no testcases" begin
    exitcode = nested_run([])
    @test exitcode == 0
end

end


text:
using DifferentialEquations
using Test
using MRIgeneralizedBloch
using MRIgeneralizedBloch: apply_hamiltonian_freeprecession!, apply_hamiltonian_gbloch!

##
max_error = 5e-2

## set parameters
ω0 = 200.0
m0s = 0.15
R1f = 0.3
R2f = 1 / 65e-3
R1s = 2
T2s = 10e-6
Rx = 30.0
TE = 3.5e-3
mfun(p, t; idxs=nothing) = typeof(idxs) <: Number ? 0.0 : zeros(30)
alg = Tsit5()

## random initialization
ms = m0s * rand()
mf = (1 - m0s) * rand()
ϑ = rand() * π / 2
φ = rand() * 2π
m0 = [mf * sin(ϑ) * cos(φ), mf * sin(ϑ) * sin(φ), mf * cos(ϑ), ms, 1]

## baseline ODE solution
FP_sol = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, m0s, R1f, R2f, Rx, R1s)), alg)

## Analytical gradients (using ApproxFun)
grad_list = [grad_m0s(), grad_R1f(), grad_R2f(), grad_Rx(), grad_R1s(), grad_T2s(), grad_ω0(), grad_B1()]

FP_sol_grad = solve(ODEProblem(apply_hamiltonian_freeprecession!, [m0; zeros(5 * (length(grad_list)),1)], (0, TE), (ω0, m0s, R1f, R2f, Rx, R1s, grad_list)), alg)

## FD derivative wrt. m0s
dm0s = 1e-9
FP_sol_dm0s = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, (m0s + dm0s), R1f, R2f, Rx, R1s)), alg)

t = 0 : 1e-5 : TE
dxf = similar(t)
dyf = similar(t)
dzf = similar(t)
dzs = similar(t)
dxf_fd = similar(t)
dyf_fd = similar(t)
dzf_fd = similar(t)
dzs_fd = similar(t)

for i = 1 : length(t)
        dxf_fd[i] = (FP_sol_dm0s(t[i])[1] - FP_sol(t[i])[1]) /dm0s
        dyf_fd[i] = (FP_sol_dm0s(t[i])[2] - FP_sol(t[i])[2]) /dm0s
        dzf_fd[i] = (FP_sol_dm0s(t[i])[3] - FP_sol(t[i])[3]) /dm0s
        dzs_fd[i] = (FP_sol_dm0s(t[i])[4] - FP_sol(t[i])[4]) /dm0s

        dxf[i] = FP_sol_grad(t[i])[6]
        dyf[i] = FP_sol_grad(t[i])[7]
        dzf[i] = FP_sol_grad(t[i])[8]
        dzs[i] = FP_sol_grad(t[i])[9]
end

@test dxf ≈ dxf_fd atol = 1e-3
@test dyf ≈ dyf_fd atol = 1e-3
@test dzf ≈ dzf_fd rtol = max_error
@test dzs ≈ dzs_fd rtol = max_error

## test derivative wrt. R1f
dR1f = 1e-9
FP_sol_dR1f = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, m0s, (R1f+ dR1f), R2f, Rx, R1s)), alg)

for i = 1: length(t)
        dxf_fd[i] = (FP_sol_dR1f(t[i])[1] - FP_sol(t[i])[1]) /dR1f
        dyf_fd[i] = (FP_sol_dR1f(t[i])[2] - FP_sol(t[i])[2]) /dR1f
        dzf_fd[i] = (FP_sol_dR1f(t[i])[3] - FP_sol(t[i])[3]) /dR1f
        dzs_fd[i] = (FP_sol_dR1f(t[i])[4] - FP_sol(t[i])[4]) /dR1f

        dxf[i] = FP_sol_grad(t[i])[11]
        dyf[i] = FP_sol_grad(t[i])[12]
        dzf[i] = FP_sol_grad(t[i])[13]
        dzs[i] = FP_sol_grad(t[i])[14]
end

@test dxf ≈ dxf_fd atol = 1e-3
@test dyf ≈ dyf_fd atol = 1e-3
@test dzf ≈ dzf_fd rtol = max_error
@test dzs ≈ dzs_fd rtol = max_error

## test derivative wrt. R2f
dR2f = 1e-6
FP_sol_dR2f = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, m0s, R1f, (R2f+dR2f), Rx, R1s)), alg)

for i = 1: length(t)
        dxf_fd[i] = (FP_sol_dR2f(t[i])[1] - FP_sol(t[i])[1]) /dR2f
        dyf_fd[i] = (FP_sol_dR2f(t[i])[2] - FP_sol(t[i])[2]) /dR2f
        dzf_fd[i] = (FP_sol_dR2f(t[i])[3] - FP_sol(t[i])[3]) /dR2f
        dzs_fd[i] = (FP_sol_dR2f(t[i])[4] - FP_sol(t[i])[4]) /dR2f

        dxf[i] = FP_sol_grad(t[i])[16]
        dyf[i] = FP_sol_grad(t[i])[17]
        dzf[i] = FP_sol_grad(t[i])[18]
        dzs[i] = FP_sol_grad(t[i])[19]
end

@test dxf ≈ dxf_fd rtol = max_error
@test dyf ≈ dyf_fd rtol = max_error
@test dzf ≈ dzf_fd atol = 1e-4
@test dzs ≈ dzs_fd atol = 1e-4

## test derivative wrt. Rx
dRx = 1e-9
FP_sol_dRx = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, m0s, R1f, R2f, (Rx+dRx), R1s)), alg)

for i = 1: length(t)
        dxf_fd[i] = (FP_sol_dRx(t[i])[1] - FP_sol(t[i])[1]) /dRx
        dyf_fd[i] = (FP_sol_dRx(t[i])[2] - FP_sol(t[i])[2]) /dRx
        dzf_fd[i] = (FP_sol_dRx(t[i])[3] - FP_sol(t[i])[3]) /dRx
        dzs_fd[i] = (FP_sol_dRx(t[i])[4] - FP_sol(t[i])[4]) /dRx

        dxf[i] = FP_sol_grad(t[i])[21]
        dyf[i] = FP_sol_grad(t[i])[22]
        dzf[i] = FP_sol_grad(t[i])[23]
        dzs[i] = FP_sol_grad(t[i])[24]
end

@test dxf ≈ dxf_fd atol = 1e-4
@test dyf ≈ dyf_fd atol = 1e-4
@test dzf ≈ dzf_fd rtol = max_error
@test dzs ≈ dzs_fd rtol = max_error

## test derivative wrt. R1s
dR1s = 1e-9
FP_sol_dR1s = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, m0s, R1f, R2f, Rx, (R1s+ dR1s))), alg)

for i = 1: length(t)
        dxf_fd[i] = (FP_sol_dR1s(t[i])[1] - FP_sol(t[i])[1]) /dR1s
        dyf_fd[i] = (FP_sol_dR1s(t[i])[2] - FP_sol(t[i])[2]) /dR1s
        dzf_fd[i] = (FP_sol_dR1s(t[i])[3] - FP_sol(t[i])[3]) /dR1s
        dzs_fd[i] = (FP_sol_dR1s(t[i])[4] - FP_sol(t[i])[4]) /dR1s

        dxf[i] = FP_sol_grad(t[i])[26]
        dyf[i] = FP_sol_grad(t[i])[27]
        dzf[i] = FP_sol_grad(t[i])[28]
        dzs[i] = FP_sol_grad(t[i])[29]
end

@test dxf ≈ dxf_fd atol = 1e-3
@test dyf ≈ dyf_fd atol = 1e-3
@test dzf ≈ dzf_fd rtol = max_error
@test dzs ≈ dzs_fd rtol = max_error

## test derivative wrt. T2s
dT2s = 1e-14
FP_sol_dT2s = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, m0s, R1f, R2f, Rx, R1s)), alg)

for i = 1: length(t)
        dxf_fd[i] = (FP_sol_dT2s(t[i])[1] - FP_sol(t[i])[1]) /dT2s
        dyf_fd[i] = (FP_sol_dT2s(t[i])[2] - FP_sol(t[i])[2]) /dT2s
        dzf_fd[i] = (FP_sol_dT2s(t[i])[3] - FP_sol(t[i])[3]) /dT2s
        dzs_fd[i] = (FP_sol_dT2s(t[i])[4] - FP_sol(t[i])[4]) /dT2s

        dxf[i] = FP_sol_grad(t[i])[31]
        dyf[i] = FP_sol_grad(t[i])[32]
        dzf[i] = FP_sol_grad(t[i])[33]
        dzs[i] = FP_sol_grad(t[i])[34]
end

@test dxf ≈ dxf_fd rtol = max_error
@test dyf ≈ dyf_fd rtol = max_error
@test dzf ≈ dzf_fd rtol = max_error
@test dzs ≈ dzs_fd rtol = max_error

## test derivative wrt. ω0
dω0 = 1
FP_sol_dω0 = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), ((ω0 + dω0), m0s, R1f, R2f, Rx, R1s)), alg)

for i = 1: length(t)
        dxf_fd[i] = (FP_sol_dω0(t[i])[1] - FP_sol(t[i])[1]) /dω0
        dyf_fd[i] = (FP_sol_dω0(t[i])[2] - FP_sol(t[i])[2]) /dω0
        dzf_fd[i] = (FP_sol_dω0(t[i])[3] - FP_sol(t[i])[3]) /dω0
        dzs_fd[i] = (FP_sol_dω0(t[i])[4] - FP_sol(t[i])[4]) /dω0

        dxf[i] = FP_sol_grad(t[i])[36]
        dyf[i] = FP_sol_grad(t[i])[37]
        dzf[i] = FP_sol_grad(t[i])[38]
        dzs[i] = FP_sol_grad(t[i])[39]
end

@test dxf ≈ dxf_fd rtol = max_error
@test dyf ≈ dyf_fd rtol = max_error
@test dzf ≈ dzf_fd atol = 1e-4
@test dzf ≈ dzf_fd atol = 1e-4

## test derivative wrt. B1
dB1 = 1e-9
FP_sol_dB1 = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, m0s, R1f, R2f, Rx, R1s)), alg)

for i = 1: length(t)
        dxf_fd[i] = (FP_sol_dB1(t[i])[1] - FP_sol(t[i])[1]) /dB1
        dyf_fd[i] = (FP_sol_dB1(t[i])[2] - FP_sol(t[i])[2]) /dB1
        dzf_fd[i] = (FP_sol_dB1(t[i])[3] - FP_sol(t[i])[3]) /dB1
        dzs_fd[i] = (FP_sol_dB1(t[i])[4] - FP_sol(t[i])[4]) /dB1

        dxf[i] = FP_sol_grad(t[i])[41]
        dyf[i] = FP_sol_grad(t[i])[42]
        dzf[i] = FP_sol_grad(t[i])[43]
        dzs[i] = FP_sol_grad(t[i])[44]
end

@test dxf ≈ dxf_fd rtol = max_error
@test dyf ≈ dyf_fd rtol = max_error
@test dzf ≈ dzf_fd rtol = max_error
@test dzs ≈ dzs_fd rtol = max_error

#############################################################
## R1a
#############################################################
R1a = 1
grad_list = [grad_R1a()]

## baseline ODE solution
FP_sol = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, m0s, R1a, R2f, Rx, R1a)), alg)

FP_sol_grad = solve(ODEProblem(apply_hamiltonian_freeprecession!, [m0; zeros(5 * (length(grad_list)),1)], (0, TE), (ω0, m0s, R1f, R2f, Rx, R1s, grad_list)), alg)

## test derivative wrt. R1a
dR1a = 1e-9
FP_sol_dR1f = solve(ODEProblem(apply_hamiltonian_freeprecession!, m0, (0, TE), (ω0, m0s, (R1a + dR1a), R2f, Rx, (R1a + dR1a))), alg)

for i = 1: length(t)
        dxf_fd[i] = (FP_sol_dR1f(t[i])[1] - FP_sol(t[i])[1]) /dR1a
        dyf_fd[i] = (FP_sol_dR1f(t[i])[2] - FP_sol(t[i])[2]) /dR1a
        dzf_fd[i] = (FP_sol_dR1f(t[i])[3] - FP_sol(t[i])[3]) /dR1a
        dzs_fd[i] = (FP_sol_dR1f(t[i])[4] - FP_sol(t[i])[4]) /dR1a

        dxf[i] = FP_sol_grad(t[i])[6]
        dyf[i] = FP_sol_grad(t[i])[7]
        dzf[i] = FP_sol_grad(t[i])[8]
        dzs[i] = FP_sol_grad(t[i])[9]
end

@test dxf ≈ dxf_fd atol = 1e-3
@test dyf ≈ dyf_fd atol = 1e-3
@test dzf ≈ dzf_fd rtol = max_error
@test dzs ≈ dzs_fd rtol = max_error

## plot the functions for debug purposes
# using Plots
# plotlyjs()
# plot(t, dxf, ticks=:native, label="xf", legend=:bottomleft)
# plot!(t, dxf_fd, label="xf_fd")
# plot!(t, dyf, ticks=:native, label="yf")
# plot!(t, dyf_fd, label="yf_fd")
# plot!(t, dzf, ticks=:native, label="zf")
# plot!(t, dzf_fd, label="zf_fd")
# plot!(t, dzs, label="zs")
# plot!(t, dzs_fd, label="zs_fd")


text:
using HorizonSideRobots
originRobot = Robot(animate = false, "24example.sit")
include("../AbstractType.jl")
include("../functions.jl")

function AverangeTemperature(robot::Robot, moveSide::HorizonSide)::Real
    averageTemperature = 0
    numOfsquares = 0
    height = 0
    stepsToBack = 0
    tempBorderSide = anti_side(moveSide)
    for _ in 1:4
        while isborder(robot, tempBorderSide)
            if (isborder(robot, tempBorderSide))
                numOfsquares += 1
                averageTemperature =+ temperature(robot)
                if (tempBorderSide == next_side_pr(next_side_pr(anti_side(moveSide))))
                    height += 1
                elseif (tempBorderSide == anti_side(moveSide))
                    stepsToBack += 1
                end
            end
            move!(robot, next_side(tempBorderSide))
        end
        move!(robot, tempBorderSide)
        tempBorderSide = next_side_pr(tempBorderSide)
    end
    for _ in 1:(height - stepsToBack)
        move!(robot, next_side_pr(moveSide))
        numOfsquares += 1
        averageTemperature =+ temperature(robot)
    end
    averageTemperature = ( averageTemperature - temperature(robot) ) / (numOfsquares - 1)
    return averageTemperature
end

function master!(robot::Robot)::Real
    maxAverageTemperature = 0
    tempAvetageTemperature = 0
    while !(isborder(robot, Nord))
        for side in [Ost, West]
            while true
                if isborder(robot, side)
                    if try_move!(robot, side)
                        tempAvetageTemperature = AverangeTemperature(robot, side)
                        if tempAvetageTemperature > maxAverageTemperature
                            maxAverageTemperature = tempAvetageTemperature
                            # можно написать очень похожую на square функцию, которая бы ещё и обвела этит прямоугольник
                        end
                    else
                        break
                    end
                else
                    move!(robot, side)
                end
            end
            if !(isborder(robot, Nord))
                move!(robot, Nord)
            end
        end
    end
    return maxAverageTemperature
end

master!(originRobot)
show!(originRobot)

max_stars_repo_name: Igetsuga/AlgorithmicsHomework

max_issues_repo_name: Igetsuga/AlgorithmicsHomework

max_forks_repo_name: Igetsuga/AlgorithmicsHomework


text:
module TableIndexes

export TableIndex

using Tables

struct TableIndex{T1, T2}
    table::T1
    colnames::Vector{Symbol}      # Column names used to construct the index
    index::Dict{T2, Vector{Int}}  # (val1, val2, ...) => [rowindex1, ...], where val_j = value of colnames[j]
end

function TableIndex(table::T1, colnames::Vector{Symbol}) where {T1}
    # Determine T2
    schema = Tables.schema(table)
    colname2type = Dict{Symbol, Any}(colname => tp for (colname, tp) in zip(schema.names, schema.types))
    types = [colname2type[colname] for colname in colnames]
    T2    = string(Tuple(types))       # "(T1, T2, ...)"
    T2    = "Tuple{$(T2[2:(end-1)])}"  # "Tuple{T1, T2, ...}"
    T2    = eval(Meta.parse(T2))

    # Construct TableIndex
    index = Dict{T2, Vector{Int}}()
    rows  = Tables.rows(table)
    row1  = Vector{Any}(undef, length(colnames))  # A row of the table
    i     = 0
    for row in rows
        i += 1
        j  = 0
        for colname in colnames
            j += 1
            row1[j] = getproperty(row, colname)
        end
        k = Tuple(row1)
        update!(index, k, i)
    end
    TableIndex{T1, T2}(table, colnames, index)
end

function update!(index, k, i)
    if !haskey(index, k)
        index[k] = [i]
    else
        push!(index[k], i)
    end
end

end


text:
__precompile__(true)

module KenyaCoVSD

using DynamicHMC,Dates,OrdinaryDiffEq,DiffEqCallbacks,JLD2,Parameters,Distributions
using NamedArrays,TransformVariables,LogDensityProblems,LinearAlgebra,BlackBoxOptim,Random
using DynamicHMC.Diagnostics,MCMCDiagnostics,MCMCChains

relative_testing_rate = vcat([0.0 for d in Date(2020,2,21):Day(1):Date(2020,3,14)],
        [t/120 for (t,d) in enumerate(Date(2020,3,15):Day(1):Date(2020,7,12))],
        [1.0 for k = 1:500])

include("transmissiondynamics.jl")
include("types.jl")
include("likelihood_functions.jl")
include("priors.jl")
include("inference.jl")

#Load defaults for the PCR and serological sensitivity after infection
@load("data/rel_sero_detection_after_infection.jld2")
@load("data/smoothedPCR_detection_after_infection.jld2")
@load("data/semi_urban_rural_detection_fits.jld2")
@load("data/rural_fits.jld2")

end # module


text:
### A Pluto.jl notebook ###
# v0.19.4

using Markdown
using InteractiveUtils

# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error).
macro bind(def, element)
    quote
        local iv = try Base.loaded_modules[Base.PkgId(Base.UUID("6e696c72-6542-2067-7265-42206c756150"), "AbstractPlutoDingetjes")].Bonds.initial_value catch; b -> missing; end
        local el = $(esc(element))
        global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : iv(el)
        el
    end
end

# ╔═╡ 61849c40-2a3b-4d53-a445-59492ad616e1
using Pkg; Pkg.activate("/Users/jj/JuliaProjects/LaserLab/")

# ╔═╡ 33dffd80-d690-11ec-0727-2552b82290be
begin
	using PlutoUI
	using CSV
	using DataFrames
	using Colors
	using Plots
	using Printf
	using Interpolations
	using QuadGK
	using Markdown
	using InteractiveUtils
	using LsqFit
	using Statistics
	using StatsBase
	using Unitful
	using UnitfulEquivalences
	using PhysicalConstants
	import Glob
end

# ╔═╡ a3e442db-eadc-470c-a0c5-a030cfdef360
import Unitful:
    nm, μm, mm, cm, m, km,
    mg, g, kg,
    ps, ns, μs, ms, s, minute, hr, d, yr, Hz, kHz, MHz, GHz,
    eV,
    μJ, mJ, J,
	μW, mW, W,
    A, N, mol, mmol, V, L, M

# ╔═╡ 2fadfcba-e14e-418a-939e-fc4894dc13b1
function ingredients(path::String)
	# this is from the Julia source code (evalfile in base/loading.jl)
	# but with the modification that it returns the module instead of the last object
	name = Symbol(basename(path))
	m = Module(name)
	Core.eval(m,
        Expr(:toplevel,
             :(eval(x) = $(Expr(:core, :eval))($name, x)),
             :(include(x) = $(Expr(:top, :include))($name, x)),
             :(include(mapexpr::Function, x) = $(Expr(:top, :include))(mapexpr, $name, x)),
             :(include($path))))
	m
end

# ╔═╡ 0674ecbc-2166-4770-8775-52de84611c2e
lfi = ingredients("../src/LaserLab.jl")

# ╔═╡ b11df77d-5d33-4393-b8d7-96c53f431f69
begin
	cmdir   = "/Users/jj/JuliaProjects/LaserLab/data/CMOS"
	odir    = "/Users/jj/JuliaProjects/LaserLab/data/CMOS"
	dfiles  = "*.csv"
	dplots  = "*.png"
	md"""
	CMOS dir = $cmdir
	"""
end

# ╔═╡ e24137d1-213b-4968-87f8-0c83858bbf69
md"""
## Single point plot
"""

# ╔═╡ d03753c6-8ed0-4c56-aa51-7e974f944274
let
	readdir(cmdir)
	dirs = lfi.LaserLab.getbolddirs(cmdir)
	md""" Select experiment : $(@bind sexp Select(dirs))"""
end

# ╔═╡ b54e7056-29a1-4a5e-928a-00a5fb1e6ad9
begin
	md""" Select point : $(@bind pnt NumberField(1:9, default=1)
)"""
end

# ╔═╡ a385e6d5-195a-474e-8a59-777b7c647d65
md"""
## Single point comparison between experiments: G2
"""

# ╔═╡ 6da703a7-20c3-4d7c-a9b6-755d5f7252ef

# ╔═╡ 6bbf2263-858c-4460-b54b-8302215dca51
md"""
## Functions
"""

# ╔═╡ 0a623b94-459e-472f-b630-6d43116420f9
function select_files(cmdir,sexp, dfiles="*.csv")
	path = joinpath(cmdir,sexp,"csv")
	readdir(path)
	xfiles = Glob.glob(dfiles, path)
	#nxfiles = string.([split(f,"/")[end] for f in xfiles])
	xfiles
end

# ╔═╡ 7c85c42f-029e-4d8c-9f17-609063f9eade
xfiles = select_files(cmdir,sexp)

# ╔═╡ d6a0e346-0e93-4f85-9f69-28d514fb2859
xfiles

# ╔═╡ e0b66717-fe68-4e18-87dd-b607ed8312cf
xfiles[pnt]

# ╔═╡ 55096351-4df1-4776-8462-ed2c37cdc997
dfpnt = DataFrame(CSV.File(xfiles[pnt]));

# ╔═╡ cd97fced-e115-4eed-8d04-bad3c6e60649
begin
	g2dirs = ["G2_BOLD_073_B1", "G2_BOLD_078_B2", "G2_BOLD_078_C2"]
	XF = [select_files(cmdir,g2d) for g2d in g2dirs]
end

# ╔═╡ ff5589b7-08ea-4911-b395-1e9b3c455c4c
function plot_spectrum_for_point(sdfp, fscale="cflt", escale="sumpes")
	plt = plot(sdfp[!, fscale], sdfp[!, escale], lw=2, label="",
		       xtickfontsize=8,ytickfontsize=8)
	scatter!(sdfp[!, fscale], sdfp[!, escale], label="")
	xlabel!("λ (nm)")
	ylabel!("pes")
	#yticks!([2e+5,4e+5,6e+5])
    xticks!([0,400,600,800])

	#xtickfontsize=18,ytickfontsize=18,xlabel="wavelength",xguidefontsize=18,yscale=:log10,ylabel="flux",yguidefontsize=18,legendfontsize=18) here

	plt
end

# ╔═╡ eb1336af-6868-4a8b-a3b7-2c7738f965e3
plot_spectrum_for_point(dfpnt)

# ╔═╡ 8ce4c1f7-fd0e-4435-9e78-38940b034890
function plot_spectrum_for_points(DF, fscale="cflt", escale="sumpes")
	plt = plot(sdfp[!, fscale], sdfp[!, escale], lw=2, label="",
		       xtickfontsize=8,ytickfontsize=8)
	scatter!(sdfp[!, fscale], sdfp[!, escale], label="")
	xlabel!("λ (nm)")
	ylabel!("pes")
	#yticks!([2e+5,4e+5,6e+5])
    xticks!([0,400,600,800])

	#xtickfontsize=18,ytickfontsize=18,xlabel="wavelength",xguidefontsize=18,yscale=:log10,ylabel="flux",yguidefontsize=18,legendfontsize=18) here

	plt
end

# ╔═╡ Cell order:
# ╠═61849c40-2a3b-4d53-a445-59492ad616e1
# ╠═33dffd80-d690-11ec-0727-2552b82290be
# ╠═a3e442db-eadc-470c-a0c5-a030cfdef360
# ╠═2fadfcba-e14e-418a-939e-fc4894dc13b1
# ╠═0674ecbc-2166-4770-8775-52de84611c2e
# ╠═b11df77d-5d33-4393-b8d7-96c53f431f69
# ╠═e24137d1-213b-4968-87f8-0c83858bbf69
# ╠═d03753c6-8ed0-4c56-aa51-7e974f944274
# ╠═7c85c42f-029e-4d8c-9f17-609063f9eade
# ╠═d6a0e346-0e93-4f85-9f69-28d514fb2859
# ╠═b54e7056-29a1-4a5e-928a-00a5fb1e6ad9
# ╠═e0b66717-fe68-4e18-87dd-b607ed8312cf
# ╠═55096351-4df1-4776-8462-ed2c37cdc997
# ╠═eb1336af-6868-4a8b-a3b7-2c7738f965e3
# ╠═a385e6d5-195a-474e-8a59-777b7c647d65
# ╠═6da703a7-20c3-4d7c-a9b6-755d5f7252ef
# ╠═cd97fced-e115-4eed-8d04-bad3c6e60649
# ╠═6bbf2263-858c-4460-b54b-8302215dca51
# ╠═0a623b94-459e-472f-b630-6d43116420f9
# ╠═ff5589b7-08ea-4911-b395-1e9b3c455c4c
# ╠═8ce4c1f7-fd0e-4435-9e78-38940b034890

max_stars_repo_name: jjgomezcadenas/LaserLab

max_issues_repo_name: jjgomezcadenas/LaserLab

max_forks_repo_name: jjgomezcadenas/LaserLab


text:
module InternalModule

function foo()
    println("InternalModule.foo()")
end

end

max_stars_repo_name: erwanlecarpentier/JuliaSandbox

max_issues_repo_name: erwanlecarpentier/JuliaSandbox

max_forks_repo_name: erwanlecarpentier/JuliaSandbox


text:
# Helper Functions for risk measures

function pnorm(p,prob,Y)
    tmp = transpose(Y)
    A = (dot(prob,tmp.^p))^(1/p)
    return A
end

function ontoSimplex(x)
    #projets given probabilities onto standard simplex to satisfy
    #sum(probabilities)= 1 and probabilities>= 0
    #works for co- and contravariant vectors
        C = x.<0                #replaces numbers under 0 with 0 in x (array)
        k = 1
        for t=1:length(x)
            if C[t] == 1
                x[t] = 0
                k = k+1
            end
        end
        tmp = sum(x)
        n = length(x)
        if tmp > 0
            x = x/tmp
        else
            #fehlt noch wird aber nicht benötigt

        end
        return x
    end

function goldenSearch(fun::Function,x1::Float64)
    #	fun: function to be minimized
    #	x1: initial value. x2: second guess (optional)
    #	example: goldenSearch(@(x)x^2, 1)
        if x1 < 0
            x2 = 1.05*x1-1e-4
        else
            x2 = 1.05*x1+1e-4
        end
        f1= fun(x1)
        f2= fun(x2)
        f3= -Inf
        phi= (sqrt(5)- 1)/2
        x3= (2+phi)*x1 - (1+phi)*x2
        while (f3 <= f1)		#finde zuerst zwei Punkte, wo dazwischen ein Minimum ist.
            x3= (2+phi)*x1 - (1+phi)*x2	#	extrapolate towards minimum by keeping golden ratio
            f3= fun(x3)
            if (f3 <= f1)		#maintain the relation f1 < f2
                x2= x1
                f2= f1
                x1= x3
                f1= f3
            else
                tmp= x1
                x1= x3
                x3= x2
                x2= tmp
                tmp= f1
                f1= f3
                f3= f2
                f2= tmp
                if (x1 > x3)
                    tmp= x1
                    x1= x3
                    x3= tmp
                end
                break
            end
        end
    #Es gelte jetzt x1 < x2 < x3 mit f1> f2 und f2< f3
    #while abs(x3- x1)> TolX*(1+ abs(x2))	%	mit golden Section search verbessere die Loesung
        #x3= (2+phi)*x1 - (1+phi)*x2
        found= false
        deltaX= Inf
        while found == false
            tmp= x3 - x1
            if tmp<deltaX  #   run to machine epsilon
                deltaX= tmp
            else
                found = true
            end
            if (x3-x2 > x2-x1)
                x = phi*x2 + (1-phi)*x3
            else
                x = (1-phi)*x1 + phi*x2
            end
            f= fun(x)
            if (f< f2)
                if (x3-x2 > x2-x1)
                    x1= x2
                    x2= x
                    f1= f2
                    f2= f
                else
                    x3= x2
                    x2= x
                    f3= f2
                    f2= f
                end
            else
                if (x3-x2 > x2-x1)
                    x3= x
                    f3= f
                else
                    x1= x
                    f1= f
                end
            end
        end

        if (f2< f1)
            tmp= f2
            f2= f1
            f1= tmp
            tmp= x2
            x2= x1
            x1= tmp
        end
        if (f3< f1)
            tmp= f3
            f3= f1
            f1= tmp
            tmp= x3
            x3= x1
            x1= tmp
        end
        #if (f3< f2)
        #    tmp= f3
        #    f3= f2
        #    f2= tmp
        #    tmp= x3
        #    x3= x2
        #    x2= tmp

        #end
        return(x1,f1)
    end

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

function checkSpectral(spec)
    if (quadgk(x -> spec(x), 0, 1, rtol=1e-8)[1] ≈ 1.0) == false
        return(false)
    end

    points = collect(range(0,stop=1,length=1000));
    values = spec(points)

    if (minimum(diff(values)) .>= 0) == true
        return(true)
    else
        return(false)
    end

end

function add_expr(op1, op2)
    expr = Expr(:call, :+, op1, op2)
    return expr
end

function math_expr(op,op1)
    expr = Expr(:call, op, op1)
    return expr
end
function math_expr(op,op1,op2)
    expr = Expr(:call, op, op1, op2)
    return expr
end