text:
module OpenImages
using HTTP
using CSV
using DataFrames
using AWSS3
using AWSCore

OID_v4 = "https://storage.googleapis.com/openimages/2018_04/"
OID_v5 = "https://storage.googleapis.com/openimages/v5/"

"""
    class_label_codes(class_labels, csv_dir=nothing)
Gets a dictionary that maps a list of OpenImages image class labels to their
corresponding image class label codes.

# Arguments
- `class_labels`: list of labels for OpenImages dataset to download.
- `csv_dir`: directory where we should look for the class annotation CSV file.
"""
function class_label_codes(class_labels, csv_dir=nothing)
    classes_csv = "class-descriptions-boxable.csv"
    url = OID_v5 * classes_csv
    response = HTTP.get(url, redirect=true)

    if response.status != 200
        println("Error!")
    end

    df_classes = CSV.File(response.body, header=false)

    if csv_dir === nothing
        csv_dir = pwd()
    end
    descriptions_csv_file_path = joinpath(csv_dir, classes_csv)

    df_classes = CSV.File(response.body, header=false)
    CSV.write(descriptions_csv_file_path, df_classes)

    labels_to_codes = Dict()
    for class_label in df_classes
        if class_label.Column2 in class_labels
            labels_to_codes[class_label.Column2] = class_label.Column1
        end
    end
    return labels_to_codes
end

"""
    get_annotations_csv(section)
Requests the annotations CSV for a split section.

# Arguments
- `section`: the relevant split section, "train", "validation", or "test".
"""
function get_annotations_csv(section)
    url = OID_v5 * section * "-annotations-bbox.csv"
    response = HTTP.get(url, redirect=true)

    if response.status != 200
        println("Error!")
    end

    return response
end

"""
    group_bounding_boxes(section, label_codes, csv_dir=nothing)
Returns a dictionary with image label as key and GroupedDataFrame as the value.

# Arguments
- `section`: the relevant split section, "train", "validation", or "test".
- `label_codes`: dictionary with class labels mapped to the image class.
- `csv_dir`: directory where we should look for the class annotation CSV file.
"""
function group_bounding_boxes(section, label_codes, csv_dir=nothing)
    if csv_dir === nothing
        response = get_annotations_csv(section)
        df_images = CSV.File(response.body)
    else
        bbox_csv_file_path = joinpath(csv_dir, section * "-annotations-bbox.csv")
        if !isfile(bbox_csv_file_path)
            response = get_annotations_csv(section)
            body = CSV.File(response.body)
            CSV.write(bbox_csv_file_path, body)
        end
        df_images = CSV.File(bbox_csv_file_path) |> DataFrame
    end

    unnecessary_columns = [
        "IsOccluded",
        "IsTruncated",
        "IsGroupOf",
        "IsDepiction",
        "IsInside",
        "Source",
        "Confidence",
    ]

    deletecols!(df_images, unnecessary_columns)

    labels_to_bounding_box_groups = Dict()
    for (labels, codes) in label_codes
        df_label_images = df_images[df_images.LabelName .== codes,:]
        deletecols!(df_label_images, ["LabelName"])
        labels_to_bounding_box_groups[labels] = groupby(df_label_images, "ImageID")
    end
    return labels_to_bounding_box_groups
end

"""
    download_single_image(image_file_path, dest_file_path)
Downloads and saves an image file from the OpenImages dataset.
"""
function download_single_image(image_file_path, dest_file_path)
    try
        config = AWSCore.aws_config(creds=nothing)
        AWSS3.s3_get_file(config, "open-images-dataset", image_file_path, dest_file_path)
    catch
        print("Error!")
    end
end

"""
    download_images_by_id(image_ids, section, image_directory)
Downloads a collection of images from OpenImages dataset.
"""
function download_images_by_id(image_ids, section, image_directory)
    for image_id in image_ids
        image_file_name = image_id * ".jpg"
        image_file_path = section * "/" * image_file_name
        dest_file_path = joinpath(image_directory, image_file_name)
        download_single_image(image_file_path, dest_file_path)
    end
end

"""
    download_images(dest_dir, class_labels, csv_dir=nothing, limit=nothing)
Downloads a dataset of images for a specified list of OpenImages image classes and returns
dictionary of the images directory for each class label.

# Arguments
- `dest_dir`: base directory under which the images and annotations will be stored.
- `class_labels`: list of labels for OpenImages dataset to download.
- `csv_dir`: directory where we should look for required CSV file (if not present files will be downloaded).
- `limit`: the maximum number of images per label to be download.
"""
function download_images(dest_dir, class_labels, csv_dir=nothing, limit=nothing)
    if csv_dir !== nothing && !isdir(csv_dir)
        mkdir(csv_dir)
    end

    labels_to_codes = class_label_codes(class_labels, csv_dir)

    class_directories = Dict()
    for (label, code) in labels_to_codes
        images_dir = joinpath(dest_dir, label, "images")
        if !isdir(images_dir)
            mkdir(images_dir)
        end

        class_directories[label] = images_dir
    end

    class_labels = [label for (label, value) in labels_to_codes]
    label_download_counts = Dict()
    for labels in class_labels
        label_download_counts[labels] = 0
    end

    for section in ("train", "validation", "test")
        label_bbox_groups = group_bounding_boxes(section, label_to_codes, csv_dir)

        for (label_index, class_label) in enumerate(class_labels)
            bbox_groups = label_bbox_groups[class_label]
            image_ids = groupvars(bbox_groups)

            if limit !== nothing
                remaining = limit - label_download_counts[class_label]
                if remaining <= 0
                    break
                elseif remaining <= length(image_ids)
                    image_ids = image_ids[1:remaining]
                end
            end
            download_images_by_id(image_ids, section, class_directories[class_label])
            label_download_counts[class_label] += length(image_ids)
        end
    end
    return class_directories
end

end # module


text:
struct SimulationStoreParams
    initial_time::Dates.DateTime
    step_resolution::Dates.Millisecond
    num_steps::Int
    # The key order is the problem execution order.
    decision_models_params::OrderedDict{Symbol, ModelStoreParams}
    emulation_model_params::OrderedDict{Symbol, ModelStoreParams}

    function SimulationStoreParams(
        initial_time::Dates.DateTime,
        step_resolution::Dates.Period,
        num_steps::Int,
        decision_models_params::OrderedDict{Symbol, ModelStoreParams},
        emulation_model_params::OrderedDict,
    )
        new(
            initial_time,
            Dates.Millisecond(step_resolution),
            num_steps,
            decision_models_params,
            emulation_model_params,
        )
    end
end

function SimulationStoreParams(initial_time, step_resolution, num_steps)
    return SimulationStoreParams(
        initial_time,
        step_resolution,
        num_steps,
        OrderedDict{Symbol, ModelStoreParams}(),
        OrderedDict{Symbol, ModelStoreParams}(),
    )
end

function SimulationStoreParams()
    return SimulationStoreParams(
        Dates.DateTime("1970-01-01T00:00:00"),
        Dates.Millisecond(0),
        0,
        OrderedDict{Symbol, ModelStoreParams}(),
        OrderedDict{Symbol, ModelStoreParams}(),
    )
end

get_initial_time(store_params::SimulationStoreParams) = store_params.initial_time

function get_decision_model_params(store_params::SimulationStoreParams, model_name::Symbol)
    return store_params.decision_models_params[model_name]
end

function get_emulation_model_params(store_params::SimulationStoreParams)
    # We currently only store one em_model dataset in the store
    @assert_op length(store_params.emulation_model_params) == 1
    return first(values(store_params.emulation_model_params))
end

ext: jl


text:
function circleofnumbers(n::Int, firstnumber::Int)::Int
    (firstnumber + n/2) % n
end

max_stars_repo_name: UlrichBerntien/Codewars-Katas

max_issues_repo_name: UlrichBerntien/Codewars-Katas

max_forks_repo_name: UlrichBerntien/Codewars-Katas


text:
__precompile__(true)
module IterativeSolvers

include("common.jl")
include("krylov.jl")
include("history.jl")

#Specialized factorizations
include("factorization.jl")

#Linear solvers
include("stationary.jl")
include("cg.jl")
include("gmres.jl")
include("chebyshev.jl")
include("idrs.jl")

#Eigensolvers
include("simple.jl")
include("lanczos.jl")

#SVD solvers
include("svdl.jl")

#Least-squares
include("lsqr.jl")
include("lsmr.jl")

#Randomized algorithms
include("rlinalg.jl")
include("rsvd.jl")
include("rsvd_fnkz.jl")

end


text:
function AUC_ROC(Ytrue, scores)
     perm = sortperm(scores)

     roc_y = Ytrue[perm]
     stack_x = cumsum(roc_y)/sum(roc_y)
     stack_y = cumsum(! roc_y)/sum(! roc_y)
     N = length(roc_y)

     auc = sum( (stack_x[2:N] - stack_x[1:N-1]) .* stack_y[2:N] )
     return auc
end


text:
module JWTs

using MbedTLS
using JSON
using Base64
using Downloads
using Random

import Base: show, isvalid
export JWT, JWK, JWKRSA, JWKSymmetric, JWKSet, issigned, isverified, isvalid, validate!, sign!, show, claims, refresh!, kid

struct JWKSymmetric
    kind::MbedTLS.MDKind
    key::Vector{UInt8}
end

struct JWKRSA
    kind::MbedTLS.MDKind
    key::Union{RSA,MbedTLS.PKContext}
end

"""
JWK represents a JWK Key (either for signing or verification).

JWK can be either a JWKRSA or JWKSymmetric. A RSA key can
represent either the public or private key.
"""
const JWK = Union{JWKRSA,JWKSymmetric}

"""
JWKSet holds a set of keys, fetched from a OpenId key URL, each key identified by a key id.

The key URL can either be of `http(s)://` or `file://` type.
"""
mutable struct JWKSet
    url::String
    keys::Dict{String,JWK}

    function JWKSet(url::String)
        new(url, Dict{String,JWK}())
    end

    function JWKSet(keyset::Vector)
        keysetdict = Dict{String,JWK}()
        refresh!(keyset, keysetdict)
        new("", keysetdict)
    end
end
function show(io::IO, jwk::JWKSet)
    print(io, "JWKSet $(length(jwk.keys)) keys")
    isempty(jwk.url) || print(io, " ($(jwk.url))")
end

"""
JWT represents a JWT payload at the minimum.

When signed, it holds the header and signature too.
The parts are stored in encoded form.
"""
mutable struct JWT
    payload::String
    header::Union{Nothing,String}
    signature::Union{Nothing,String}
    verified::Bool
    valid::Union{Nothing,Bool}

    function JWT(; jwt::Union{Nothing,String}=nothing, payload::Union{Nothing,Dict{String,Any},String}=nothing)
        if jwt !== nothing
            @assert payload === nothing
            parts = split(jwt, ".")
            if length(parts) == 3
                new(parts[2], parts[1], parts[3], false, nothing)
            else
                new("", nothing, nothing, true, false)
            end
        else
            @assert payload !== nothing
            new(isa(payload, String) ? payload : urlenc(base64encode(JSON.json(payload))), nothing, nothing, false, nothing)
        end
    end
end

decodepart(encoded::String) = JSON.parse(String(base64decode(urldec(encoded))))

claims(jwt::JWT) = decodepart(jwt.payload)
issigned(jwt::JWT) = (nothing !== jwt.signature) && (nothing !== jwt.header)
isverified(jwt::JWT) = jwt.verified
isvalid(jwt::JWT) = jwt.valid
function kid(jwt::JWT)
    @assert issigned(jwt)
    decodepart(jwt.header)["kid"]
end

show(io::IO, jwt::JWT) = print(io, issigned(jwt) ? join([jwt.header, jwt.payload, jwt.signature], '.') : jwt.payload)

validate!(jwt::JWT, keyset::JWKSet) = validate!(jwt, keyset, kid(jwt))
function validate!(jwt::JWT, keyset::JWKSet, kid::String)
    isverified(jwt) && (return isvalid(jwt))
    (kid in keys(keyset.keys)) || refresh!(keyset)
    validate!(jwt, keyset.keys[kid])
end
function validate!(jwt::JWT, key::T) where {T <: JWK}
    isverified(jwt) && (return isvalid(jwt))
    @assert issigned(jwt)

    data = jwt.header * "." * jwt.payload
    sigbytes = base64decode(urldec(jwt.signature))

    jwt.verified = true

    jwt.valid = if T <: JWKRSA
        try
            MbedTLS.verify(key.key, key.kind, MbedTLS.digest(key.kind, data), sigbytes) == 0
        catch
            false
        end
    else
        MbedTLS.digest(key.kind, data, key.key) == sigbytes
    end
end

function sign!(jwt::JWT, keyset::JWKSet, kid::String)
    issigned(jwt) && return
    (kid in keys(keyset.keys)) || refresh!(keyset)
    sign!(jwt::JWT, keyset.keys[kid], kid)
end
function sign!(jwt::JWT, key::T, kid::String="") where {T <: JWK}
    issigned(jwt) && return

    if T <: JWKRSA
        if key.kind === MbedTLS.MD_SHA256
            alg = "RS256"
        elseif key.kind === MbedTLS.MD_SHA384
            alg = "RS384"
        elseif key.kind === MbedTLS.MD_SHA
            alg = "RS512"
        else
            error("unsupported key algorithm")
        end
    else
        if key.kind === MbedTLS.MD_SHA256
            alg = "HS256"
        elseif key.kind === MbedTLS.MD_SHA384
            alg = "HS384"
        elseif key.kind === MbedTLS.MD_SHA
            alg = "HS512"
        else
            error("unsupported key algorithm")
        end
    end
    alg = (T <: JWKRSA) ? "RS256" : "HS256"
    header_dict = Dict{String,String}("alg"=>alg, "typ"=>"JWT")
    isempty(kid) || (header_dict["kid"] = kid)
    header = urlenc(base64encode(JSON.json(header_dict)))

    data = header * "." * jwt.payload
    sigbytes = (T <: JWKRSA) ?  MbedTLS.sign(key.key, key.kind, MbedTLS.digest(key.kind, data), MersenneTwister(0)) : MbedTLS.digest(key.kind, data, key.key)
    signature = urlenc(base64encode(sigbytes))

    jwt.header = header
    jwt.signature = signature
    jwt.verified = true
    jwt.valid = true
    nothing
end

function refresh!(keyset::JWKSet, keyseturl::String)
    keyset.url = keyseturl
    refresh!(keyset)
end

function refresh!(keyset::JWKSet)
    if !isempty(keyset.url)
        keys = Dict{String,JWK}()
        refresh!(keyset.url, keys)
        keyset.keys = keys
    end
    nothing
end

function refresh!(keyseturl::String, keysetdict::Dict{String,JWK})
    if startswith(keyseturl, "file://")
        jstr = readchomp(keyseturl[8:end])
    else
        output = PipeBuffer()
        Downloads.request(keyseturl; method="GET", output=output)
        jstr = String(take!(output))
    end
    keys = JSON.parse(jstr)["keys"]
    refresh!(keys, keysetdict)
end

function refresh!(keys::Vector, keysetdict::Dict{String,JWK})
    for key in keys
        kid = key["kid"]
        kty = key["kty"]
        alg = key["alg"]

        # ref: https://tools.ietf.org/html/rfc7518
        try
            if kty == "RSA"
                n = base64decode(urldec(key["n"]))
                e = base64decode(urldec(key["e"]))
                if alg == "RS256"
                    keysetdict[kid] = JWKRSA(MbedTLS.MD_SHA256, pubkey(n, e, MbedTLS.MD_SHA256))
                elseif alg == "RS384"
                    keysetdict[kid] = JWKRSA(MbedTLS.MD_SHA384, pubkey(n, e, MbedTLS.MD_SHA384))
                elseif alg == "RS512"
                    keysetdict[kid] = JWKRSA(MbedTLS.MD_SHA, pubkey(n, e, MbedTLS.MD_SHA))
                else
                    @warn("key alg $(key["alg"]) not supported yet, skipping key $kid")
                    continue
                end
            elseif kty == "oct"
                k = base64decode(urldec(key["k"]))
                if alg == "HS256"
                    keysetdict[kid] = JWKSymmetric(MbedTLS.MD_SHA256, k)
                elseif alg == "HS384"
                    keysetdict[kid] = JWKSymmetric(MbedTLS.MD_SHA384, k)
                elseif alg == "HS512"
                    keysetdict[kid] = JWKSymmetric(MbedTLS.MD_SHA, k)
                else
                    @warn("key alg $(key["alg"]) not supported yet, skipping key $kid")
                    continue
                end
            else
                @warn("key type $(key["kty"]) not supported yet, skipping key $kid")
                continue
            end
        catch ex
            @warn("exception $ex trying to decode, skipping key $kid")
        end
    end
    nothing
end

function pubkey(bytesn, bytese, halg)
    n = parse(BigInt, bytes2hex(bytesn); base=16)
    e = parse(BigInt, bytes2hex(bytese); base=16)

    R = RSA(MbedTLS.MBEDTLS_RSA_PKCS_V15, halg)
    MbedTLS.pubkey_from_vals!(R, e, n)
end

function urldec(bs)
    bs = replace(bs, "-"=>"+")
    bs = replace(bs, "_"=>"/")
    padb64(bs)
end

function urlenc(bs)
    bs = replace(bs, "+"=>"-")
    bs = replace(bs, "/"=>"_")
    bs = replace(bs, "="=>"")
    bs
end

function padb64(bs)
    surplus = length(bs) % 4
    if surplus > 0
        bs = bs * "="^(4 - surplus)
    end
    bs
end

end # module JWTs


text:
# functions to read data
# originally in functions.jl
# Claudia March 2015
# based on Julia 0.4
# Claudia October 2015

# ----- read data --------

# function to write a csv table from the expCF of an
# array of quartets
# warning: does not check if the expCF have been calculated
function writeExpCF(quartets::Array{Quartet,1})
    df = DataFrames.DataFrame(t1=String[],t2=String[],t3=String[],t4=String[],
                              CF12_34=Float64[],CF13_24=Float64[],CF14_23=Float64[])
    for q in quartets
        length(q.taxon) == 4 || error("quartet $(q.number) does not have 4 taxa")
        length(q.qnet.expCF) == 3 || error("quartet $(q.number) does have qnet with 3 expCF")
        push!(df, [q.taxon[1],q.taxon[2],q.taxon[3],q.taxon[4],q.qnet.expCF[1],q.qnet.expCF[2],q.qnet.expCF[3]])
    end
    return df
end

writeExpCF(d::DataCF) = writeExpCF(d.quartet)

# function to write a csv table from the obsCF of an
# array of quartets
function writeObsCF(quartets::Array{Quartet,1})
    df = DataFrames.DataFrame(t1=String[],t2=String[],t3=String[],t4=String[],
                              CF12_34=Float64[],CF13_24=Float64[],CF14_23=Float64[],ngenes=Int[])
    for q in quartets
        length(q.taxon) == 4 || error("quartet $(q.number) does not have 4 taxa")
        length(q.obsCF) == 3 || error("quartet $(q.number) does have qnet with 3 expCF")
        push!(df, [q.taxon[1],q.taxon[2],q.taxon[3],q.taxon[4],q.obsCF[1],q.obsCF[2],q.obsCF[3],q.ngenes])
    end
    return df
end

writeObsCF(d::DataCF) = writeObsCF(d.quartet)

# function that takes a dataframe and creates a DataCF object
# has ! because it can modify the dataframe inside
"""
    readTableCF(file)
    readTableCF(data frame)
    readTableCF!(data frame)

Read a file or DataFrame object containing a table of concordance factors (CF),
with one row per 4-taxon set. The first 4 columns are assumed to give the labels
of the 4 taxa in each set (tx1, tx2, tx3, tx4).
Columns containing the CFs are assumed to be named
'CF12_34', 'CF13_24' and 'CF14_23', or else are assumed to be columns 5,6,7.
If present, a column named 'ngenes' will be used to get the number of loci
used to estimate the CFs for each 4-taxon set.

Optional arguments:

- summaryfile: if specified, a summary file will be created with that name.
- delim (for the first form only): to specify how columns are delimited,
  with single quotes: delim=';'. Default is a `csv` file, i.e. `delim=','`.

The last version modifies the input data frame, if species are represented by multiple alleles
for instance (see `readTableCF!`).
"""
# warning: file AbstractString bc it can be read as UTF8String
function readTableCF(file::AbstractString; delim=','::Char, summaryfile=""::AbstractString)
    df = CSV.read(file, delim=delim)
    readTableCF!(df, summaryfile=summaryfile)
end

function readTableCF(df0::DataFrames.DataFrame; summaryfile=""::AbstractString)
    df = deepcopy(df0)
    readTableCF!(df, summaryfile=summaryfile)
end

function readTableCF!(df::DataFrames.DataFrame; summaryfile=""::AbstractString)
    DEBUG && println("assume the numbers for the taxon read from the observed CF table match the numbers given to the taxon when creating the object network")
    alternativecolnames = [ # obsCF12 is as exported by fittedQuartetCF()
        [:CF12_34, Symbol("CF12.34"), :obsCF12],
        [:CF13_24, Symbol("CF13.24"), :obsCF13],
        [:CF14_23, Symbol("CF14.23"), :obsCF14]
    ]
    obsCFcol = [findfirst(x-> x ∈ alternativecolnames[1], DataFrames.names(df)),
                findfirst(x-> x ∈ alternativecolnames[2], DataFrames.names(df)),
                findfirst(x-> x ∈ alternativecolnames[3], DataFrames.names(df))]
    ngenecol =  findfirst(DataFrames.names(df), :ngenes)
    withngenes = ngenecol>0
    if findfirst(obsCFcol, 0) > 0 # one or more col names for CFs were not found
        size(df,2) == (withngenes ? 8 : 7) ||
          warn("""Column names for quartet concordance factors (CFs) were not recognized.
          Was expecting CF12_34, CF13_24 and CF14_23 for the columns with CF values,
          or CF12.34 or obsCF12, etc.
          Will assume that the first 4 columns give the taxon names, and that columns 5-7 give the CFs.""")
        obsCFcol = [5,6,7] # assuming CFs are in columns 5,6,7, with colname mismatch
    end
    minimum(obsCFcol) > 4 ||
        error("CFs found in columns $obsCFcol, but taxon labels expected in columns 1-4")
    # fixit: what about columns giving the taxon names: always assumed to be columns 1-4? No warning if not?
    columns = [[1,2,3,4]; obsCFcol]
    if withngenes  push!(columns, ngenecol)  end

    d = readTableCF!(df, columns)

    if withngenes && d.numTrees == -1
        m1 = minimum([q.ngenes for q in d.quartet])
        m2 = maximum([q.ngenes for q in d.quartet])
        if m1<m2 print("between $m1 and ") end
        println("$m2 gene trees per 4-taxon set")
        # other info printed by show() on a DataCF object: num quartets and num gene trees
    end
    if(summaryfile != "")
        descData(d,summaryfile)
    end
    return d
end

# see docstring below, for readTableCF!
# takes in df and 7 or 8 column numbers (4 labels + 3 CFs + ngenes possibly)
function readTableCF!(df::DataFrames.DataFrame, co::Vector{Int})
    withngenes = (length(co)==8) # true if column :ngenes exists, false ow
    repSpecies = cleanAlleleDF!(df,co) # removes uninformative rows from df (not df0)
    # fixit: cleanAlleleDF! is time consuming but many times not needed
    # add option to skip it, if the user knows that each tip appears once only?
    if !isempty(repSpecies)
        df = mergeRows(df,co)   # warning: this 'df' is *not* changed externally
    end                         # we cannot move to mapAllelesCFtable because we need repSpecies in here
    quartets = Quartet[]
    for i in 1:size(df,1)
        push!(quartets,Quartet(i,string(df[i,co[1]]),string(df[i,co[2]]),string(df[i,co[3]]),string(df[i,co[4]]),
                               [df[i,co[5]],df[i,co[6]],df[i,co[7]]]))
        if(withngenes)
            quartets[end].ngenes = df[i,co[8]]
        end
    end
    d = DataCF(quartets)
    if(!isempty(repSpecies))
        d.repSpecies = repSpecies
    end
    return d  # return d, df ## to save memory & gc with readTableCF! for bootstrapping?
end

"""
    readTableCF!(data frame, columns)

Read in quartet CFs from data frame, assuming information is in columns numbered `columns`,
of length **7 or 8**: 4 taxon labels then 3 CFs then ngenes possibly.

If some species appears more than once in the same 4-taxon set (e.g. t1,t1,t2,t3),
then the data frame is modified to remove rows (4-taxon sets) that are uninformative about
between-species relationships. This situation may occur if multiple individuals are sampled from
the same species. A 4-taxon set is uninformative (and its row is removed)
if one taxon is repeated 3 or 4 times (like t1,t1,t1,t1 or t1,t2,t2,t2).
The list of species appearing twice in some 4-taxon sets is stored in the output DataCF object.
For these species, the length of their external edge is identifiable (in coalescent units).
If multiple rows correspond to the same 4-taxon set, these rows are merged and their CF values
(and number of genes) are averaged.

    readTableCF!(DataCF, data frame, columns)

Modify the `.quartet.obsCF` values in the `DataCF` object with those read from the data frame
in columns numbered `columns`.
`columns` should have **3** columns numbers for the 3 CFs in this order: 12_34, 13_24 and 14_23.
Assume the same 4-taxon sets in `DataCF` and in the data frame, and in the same order,
but this assumption is *not checked* (for speed, e.g. during bootstrapping).
"""
# WARNING: assumes df has a single row per 4-taxon sets (multiple ind already merged)
#                  dcf created earlier from df.
function readTableCF!(datcf::DataCF, df::DataFrame, cols::Vector{Int})
    for i in 1:size(df,1)
        for j in 1:3
            datcf.quartet[i].obsCF[j] = df[i,cols[j]]
        end
    end
end

# ---------------- read input gene trees and calculate obsCF ----------------------

# function to read a file and create one object per line read
# (each line starting with "(" will be considered a topology)
# the file can have extra lines that are ignored
# returns an array of HybridNetwork objects (that can be trees)
"""
`readInputTrees(file)`

Read a text file with a list of trees/networks in parenthetical format
(one tree per line) and transform them like [`readTopologyLevel1`](@ref)
does: to be unrooted, with resolved polytomies, missing branch lengths
set to 1.0, etc. See [`readMultiTopology`](@ref) to read multiple
trees or networks with no modification.

Output: array of HybridNetwork objects.
"""
function readInputTrees(file::AbstractString)
    try
        s = open(file)
    catch
        error("Could not find or open $(file) file");
    end
    vnet = HybridNetwork[];
    s = open(file)
    numl = 1
    for line in eachline(s)
        line = strip(line) # remove spaces
        DEBUG && println("$(line)")
        c = isempty(line) ? "" : line[1]
        if(c == '(')
           try
               push!(vnet, readTopologyUpdate(line,false))
           catch(err)
               error("could not read tree in line $(numl). The error is $(err)")
           end
        end
        numl += 1
    end
    close(s)
    return vnet # consistent output type: HybridNetwork vector. might be of length 0.
end

# function to list all quartets for a set of taxa names
# return a vector of quartet objects, and if writeFile=true, writes a file
# warning: taxon has to be vector of String, vector of AbstractString do not work
function allQuartets(taxon::Union{Vector{String},Vector{Int}}, writeFile::Bool)
    quartets = combinations(taxon,4)
    vquartet = Quartet[];
    if(writeFile)
        #allName = "allQuartets$(string(integer(time()/1000))).txt"
        allName = "allQuartets.txt"
        f = open(allName,"w")
    end
    i = 1
    for q in quartets
        if(writeFile)
            write(f,"$(q[1]),$(q[2]),$(q[3]),$(q[4])\n")
        end
        push!(vquartet,Quartet(i,string(q[1]),string(q[2]),string(q[3]),chomp(string(q[4])),[1.0,0.0,0.0]))
        i += 1 # overflow error if # quartets > typemax(Int), i.e. if 121,978+ taxa with Int64, 478+ taxa with Int32
    end
    if(writeFile)
        close(f)
    end
    return vquartet
end

allQuartets(numTaxa::Integer, writeFile::Bool) = allQuartets(1:numTaxa, writeFile)

# function to list num randomly selected quartets for a vector of all quartets
# return a vector of Quartet randomly chosen (and a file if writeFile=true)
function randQuartets(allquartets::Vector{Quartet},num::Integer, writeFile::Bool)
    randquartets = Quartet[]
    n = length(allquartets)
    if(num == 0)
        num = integer(floor(0.1*n))
    end
    num <= n || error("you cannot choose a sample of $(num) quartets when there are $(n) in total")
    indx = [ones(Int,num); zeros(Int,n-num)]
    indx = indx[sortperm(randn(n))]
    if(writeFile)
        #randName = "rand$(numQ)Quartets$(string(integer(time()/1000))).txt"
        randName = "rand$(numQ)Quartets.txt"
        println("list of randomly selected quartets in file $(randName)")
        out = open(randName,"w")
    end
    for i in 1:n
        if(indx[i] == 1)
            if(writeFile)
                write(out,"$(lines[i])")
            end
            push!(randquartets,allquartets[i])
        end
    end
    if(writeFile)
        close(out)
    end
    return randquartets
end

# function that will not use randQuartets(list of quartets,...)
# this function uses whichQuartet to avoid making the list of all quartets
# fixit: i think we should write always the file, but not sure
function randQuartets(taxon::Union{Vector{String},Vector{Int}},num::Integer, writeFile::Bool)
    randquartets = Quartet[]
    n = length(taxon)
    ntotal = binom(n,4)
    num <= ntotal || error("you cannot choose a sample of $(num) quartets when there are $(ntotal) in total")
    # indx = [rep(1,num);rep(0,ntotal-num)] # requires much more memory than necessary:
    # indx = indx[sortperm(randn(ntotal))]  # several arrays of size ntotal !!
    # rq = find(indx .== 1)
    rq = sample(1:ntotal, num, replace=false, ordered=true)
    randName = "rand$(num)Quartets.txt"
    println("list of randomly selected quartets in file $(randName)")
    out = open(randName,"w")
    i = 1
    for q in rq
        qind = whichQuartet(n,q) # vector of int
        quartet = createQuartet(taxon,qind,i)
        write(out,"$(quartet.taxon[1]), $(quartet.taxon[2]), $(quartet.taxon[3]), $(quartet.taxon[4])\n")
        push!(randquartets,quartet)
        i += 1
    end
    close(out)
    return randquartets
end

randQuartets(numTaxa::Integer,num::Integer, writeFile::Bool) = randQuartets(1:numTaxa,num, writeFile)

# function to read list of quartets from a file
# and create Quartet type objects
function readListQuartets(file::AbstractString)
    try
        f = open(file)
    catch
        error("could not open file $(file)")
    end
    f = open(file)
    quartets = Quartet[];
    i = 1
    for line in eachline(f)
        l = split(line,",")
        length(l) == 4 || error("quartet with $(length(l)) elements, should be 4: $(line)")
        push!(quartets,Quartet(i,string(l[1]),string(l[2]),string(l[3]),chomp(string(l[4])),[1.0,0.0,0.0]))
        i += 1
    end
    return quartets
end

"""
    sameTaxa(Quartet, HybridNetwork)

Return `true` if all taxa in the quartet are represented in the network,
`false` if one or more taxa in the quartet does not appear in the network.

warning: the name can cause confusion. A more appropriate name might be
"in", or "taxain", or "taxonsubset", or etc.
"""
function sameTaxa(q::Quartet, t::HybridNetwork)
    for name in q.taxon
        in(name,t.names) || return false
    end
    return true
end

"""
    taxadiff(Vector{Quartet}, network; multiplealleles=true)
    taxadiff(DataCF, network; multiplealleles=true)

Return 2 vectors:

- taxa in at least 1 of the quartets but not in the network, and
- taxa in the network but in none of the quartets.

When `multiplealleles` is true, the taxon names that end with "__2"
are ignored in the quartets: they are not expected to appear in the
networks that users give as input, or get as output.
"""
function taxadiff(quartets::Vector{Quartet}, t::HybridNetwork;
                  multiplealleles=true::Bool)
    tq = tipLabels(quartets)
    secondallele = ismatch.(r"__2$", tq)
    for i in length(secondallele):-1:1
        secondallele[i] || continue
        basetax = match(r"(.*)__2$", tq[i]).captures[1]
        # if tq[i] = "mouse__2" for instance, then basetax = "mouse"
        if findfirst(tq, basetax) > 0 # some other taxon is "mouse"
            deleteat!(tq, i) # delete "mouse__2" from tq IF "mouse" is present
        end
    end
    tn = tipLabels(t)
    return (setdiff(tq,tn), setdiff(tn,tq))
end

taxadiff(d::DataCF, t::HybridNetwork; multiplealleles=true::Bool) =
    taxadiff(d.quartet, t; multiplealleles=multiplealleles)

# function to extract the union of taxa of list of gene trees
function unionTaxa(trees::Vector{HybridNetwork})
    taxa = union([tipLabels(t) for t in trees]...) # vector, order maintained
    return sortUnionTaxa!(taxa)
end

"""
    sortUnionTaxa!(taxa)

Take a vector of strings `taxa`, sort it numerically if
elements can be parsed as an integer, alphabetically otherwise.
"""
function sortUnionTaxa!(taxa)
    integers = true
    try
        parse.(Int,taxa)
    catch
        integers = false
    end
    if integers
        sort!(taxa, by=x->parse(Int,x))
    else
        sort!(taxa)
    end
    return taxa
end

# function to extract the union of taxa of list of quartets
function unionTaxa(quartets::Vector{Quartet})
    taxa = union([q.taxon for q in quartets]...) # vector, order maintained
    return sortUnionTaxa!(taxa)
end

unionTaxaTree(file::AbstractString) = unionTaxa(readInputTrees(file))

tipLabels(t::Vector{HybridNetwork}) = unionTaxa(t)
tipLabels(q::Vector{Quartet}) = unionTaxa(q)
tipLabels(d::DataCF) = unionTaxa(d.quartet)

"""
    calculateObsCFAll!(DataCF, taxa::Union{Vector{String}, Vector{Int}})

update the .quartet[i].obsCF values of the DataCF object, based on its .tree vector.

    calculateObsCFAll!(vector of quartets, vector of trees, taxa)

update the .obsCF values of the quartets, based on the trees, and returns a new DataCF object
with these updated quartets and trees.

    calculateObsCFAll_noDataCF!(vector of quartets, vector of trees, taxa)

update the .obsCF values of the quartets based on the trees, but returns nothing.
"""
# function to calculate the obsCF from a file with a set of gene trees
# returns a DataCF object (*no* longer writes a csv table with the obsCF)
# warning: it needs trees (not networks) as input
function calculateObsCFAll!(dat::DataCF, taxa::Union{Vector{String}, Vector{Int}})
    calculateObsCFAll_noDataCF!(dat.quartet, dat.tree, taxa)
end

function calculateObsCFAll!(quartets::Vector{Quartet}, trees::Vector{HybridNetwork}, taxa::Union{Vector{String}, Vector{Int}})
    calculateObsCFAll_noDataCF!(quartets, trees, taxa)
    d = DataCF(quartets,trees)
    return d
end

function calculateObsCFAll_noDataCF!(quartets::Vector{Quartet}, trees::Vector{HybridNetwork}, taxa::Union{Vector{String}, Vector{Int}})
    println("calculating obsCF from $(length(trees)) gene trees and for $(length(quartets)) quartets")
    index = 1
    totalq = length(quartets)
    println("Reading in quartets...")
    r = round(1/totalq,2)
    numq = (r > 0.02 ? totalq : 50)
    print("0+")
    for i in 1:numq
        print("-")
    end
    print("+100%")
    println("  ")
    print("  ")
    for q in quartets
        if(round(index/totalq,2)>0.02)
            print("*")
            index = 1
        end
        suma = 0
        sum12 = 0
        sum13 = 0
        sum14 = 0
        for t in trees
            isTree(t) || error("gene tree found in file that is a network $(writeTopology(t))")
            if sameTaxa(q,t)
                M = tree2Matrix(t,taxa) #fixit: way to reuse M? length(t.edge) will be different across trees
                res = extractQuartetTree(q,M,taxa)
                DEBUG && println("res is $(res)")
                if(res == 1)
                    sum12 += 1
                elseif(res == 2)
                    sum13 += 1
                elseif(res == 3)
                    sum14 += 1
                end
                suma += (res == 0) ? 0 : 1
            end
        end
        q.obsCF = [sum12/suma, sum13/suma, sum14/suma]
        q.ngenes = suma
        index += 1
    end
    println("  ")
    return nothing
end

# function to read input list of gene trees/quartets and calculates obsCF
# as opposed to readTableCF that read the table of obsCF directly
# input: treefile (with gene trees), quartetfile (with list of quartets),
# whichQ (:add/:rand to decide if all or random sample of quartets, default all)
# numQ: number of quartets in random sample
# writetab = true to write the table of obsCF as file with name filename
# does it by default
# writeFile=true writes file with sampled quartets, default false
function readInputData(treefile::AbstractString, quartetfile::AbstractString, whichQ::Symbol, numQ::Integer, writetab::Bool, filename::AbstractString, writeFile::Bool, writeSummary::Bool)
    if writetab
        if(filename == "none")
            filename = "tableCF.txt" # "tableCF$(string(integer(time()/1000))).txt"
        end
        if (isfile(filename) && filesize(filename) > 0)
           error("""file $(filename) already exists and is non-empty. Cannot risk to erase data.
                    Choose a different CFfile name, use writeTab=false, or read the existing file
                    with readTableCF(\"$(filename)\")""")
        end
    end
    println("read input trees from file $(treefile)\nand quartetfile $(quartetfile)")
    trees = readInputTrees(treefile)
    readInputData(trees, quartetfile, whichQ, numQ, writetab, filename, writeFile, writeSummary)
end

readInputData(treefile::AbstractString, quartetfile::AbstractString, whichQ::Symbol, numQ::Integer, writetab::Bool) = readInputData(treefile, quartetfile, whichQ, numQ, writetab, "none", false, true)
readInputData(treefile::AbstractString, quartetfile::AbstractString, whichQ::Symbol, numQ::Integer) = readInputData(treefile, quartetfile, whichQ, numQ, true, "none", false, true)
##readInputData(treefile::AbstractString, quartetfile::AbstractString) = readInputData(treefile, quartetfile, :all, 0, true, "none", false, true)
readInputData(treefile::AbstractString, quartetfile::AbstractString, writetab::Bool, filename::AbstractString) = readInputData(treefile, quartetfile, :all, 0, writetab, filename, false, true)

# function to read input list of gene trees/quartets and calculates obsCF
# as opposed to readTableCF that read the table of obsCF directly
# input: trees Vector of HybridNetwork, quartetfile (with list of quartets),
# whichQ (:add/:rand to decide if all or random sample of quartets, default all)
# numQ: number of quartets in random sample
# writetab = true to write the table of obsCF as file with name filename
# does it by default
# writeFile=true writes file with sampled quartets, default false
function readInputData(trees::Vector{HybridNetwork}, quartetfile::AbstractString, whichQ::Symbol, numQ::Integer, writetab::Bool, filename::AbstractString, writeFile::Bool, writeSummary::Bool)
    if(whichQ == :all)
        numQ == 0 || warn("set numQ=$(numQ) but whichQ is not rand, so all quartets will be used and numQ will be ignored. If you want a specific number of 4-taxon subsets not random, you can input with the quartetfile option")
        println("will use all quartets in file $(quartetfile)")
        quartets = readListQuartets(quartetfile)
    elseif(whichQ == :rand)
        if(numQ == 0)
            warn("not specified numQ but whichQ=rand, so 10% of quartets will be sampled") #handled inside randQuartets
        else
            println("will take a random sample of $(numQ) 4-taxon sets from file $(quartetfile)")
        end
        allquartets = readListQuartets(quartetfile)
        quartets = randQuartets(allquartets,numQ,writeFile)
    else
        error("unknown symbol for whichQ $(whichQ), should be either all or rand")
    end
    d = calculateObsCFAll!(quartets,trees, unionTaxa(trees))
    if(writetab)
        if(filename == "none")
            filename = "tableCF.txt" # "tableCF$(string(integer(time()/1000))).txt"
        end
        if (isfile(filename) && filesize(filename) > 0)
           error("""file $(filename) already exists and is non-empty. Cannot risk to erase data.
                    Choose a different CFfile name, use writeTab=false, or read the existing file
                    with readTableCF(\"$(filename)\")""")
        end
        println("\ntable of obsCF printed to file $(filename)")
        df = writeObsCF(d)
        CSV.write(filename,df)
    end
    #descData(d,"summaryTreesQuartets$(string(integer(time()/1000))).txt")
    writeSummary && descData(d,"summaryTreesQuartets.txt")
    return d
end

# function to read input list of gene trees, and not the list of quartets
# so it creates the list of quartets inside and calculates obsCF
# as opposed to readTableCF that read the table of obsCF directly
# input: treefile (with gene trees), whichQ (:add/:rand to decide if all or random sample of quartets, default all)
# numQ: number of quartets in random sample
# taxa: list of taxa, if not given, all taxa in gene trees used
# writetab = true to write the table of obsCF as file with name filename
# does it by default
# writeFile= true, writes intermediate files with the quartets info (default false)
function readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, taxa::Union{Vector{String}, Vector{Int}}, writetab::Bool, filename::AbstractString, writeFile::Bool, writeSummary::Bool)
    if writetab
        if(filename == "none")
            filename = "tableCF.txt" # "tableCF$(string(integer(time()/1000))).txt"
        end
        if (isfile(filename) && filesize(filename) > 0)
           error("""file $(filename) already exists and is non-empty. Cannot risk to erase data.
                    Choose a different CFfile name, use writeTab=false, or read the existing file
                    with readTableCF(\"$(filename)\")""")
        end
    end
    println("read input trees from file $(treefile). no quartet file given.")
    trees = readInputTrees(treefile)
    readInputData(trees, whichQ, numQ, taxa, writetab, filename, writeFile, writeSummary)
end

readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, taxa::Union{Vector{String}, Vector{Int}}, writetab::Bool) = readInputData(treefile, whichQ, numQ, taxa, writetab, "none", false, true)
readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, taxa::Union{Vector{String}, Vector{Int}}) = readInputData(treefile, whichQ, numQ, taxa, true, "none",false, true)
readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, writetab::Bool, filename::AbstractString) = readInputData(treefile, whichQ, numQ, unionTaxaTree(treefile), writetab, filename,false, true)
readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, writetab::Bool) = readInputData(treefile, whichQ, numQ, unionTaxaTree(treefile), writetab, "none",false, true)
readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer) = readInputData(treefile, whichQ, numQ, unionTaxaTree(treefile), true, "none",false, true)
readInputData(treefile::AbstractString) = readInputData(treefile, :all, 0, unionTaxaTree(treefile), true, "none",false, true)
readInputData(treefile::AbstractString,taxa::Union{Vector{String}, Vector{Int}}) = readInputData(treefile, :all, 0, taxa, true, "none",false, true)
# above: the use of unionTaxaTree to set the taxon set
#        is not good: need to read the tree file twice: get the taxa, then get the trees
#        this inefficiency was fixed in readTrees2CF

# function to read input vector of HybridNetworks, and not the list of quartets
# so it creates the list of quartets inside and calculates obsCF
# as opposed to readTableCF that read the table of obsCF directly
# input: trees (Vector of HybridNetwork), whichQ (:add/:rand to decide if all or random sample of quartets, default all)
# numQ: number of quartets in random sample
# taxa: list of taxa, if not given, all taxa in gene trees used
# writetab = true to write the table of obsCF as file with name filename
# does it by default
# writeFile= true, writes intermediate files with the quartets info (default false)
function readInputData(trees::Vector{HybridNetwork}, whichQ::Symbol, numQ::Integer, taxa::Union{Vector{String}, Vector{Int}}, writetab::Bool, filename::AbstractString, writeFile::Bool, writeSummary::Bool)
    if(whichQ == :all)
        numQ == 0 || warn("set numQ=$(numQ) but whichQ=all, so all quartets will be used and numQ will be ignored. If you want a specific number of 4-taxon subsets not random, you can input with the quartetfile option")
        quartets = allQuartets(taxa,writeFile)
        println("will use all quartets on $(length(taxa)) taxa")
    elseif(whichQ == :rand)
        if(numQ == 0)
            warn("not specified numQ with whichQ=rand, so 10% of quartets will be sampled") #handled inside randQuartets
        else
            println("will use a random sample of $(numQ) 4-taxon sets ($(round((100*numQ)/binomial(length(taxa),4),2)) percent) on $(length(taxa)) taxa")
        end
        quartets = randQuartets(taxa,numQ, writeFile)
    else
        error("unknown symbol for whichQ $(whichQ), should be either all or rand")
    end
    d = calculateObsCFAll!(quartets,trees,taxa)
    if writetab
        if(filename == "none")
            filename = "tableCF.txt"
        end
        println("table of obsCF printed to file $(filename)")
        df = writeObsCF(d)
        CSV.write(filename,df)
    end
    #descData(d,"summaryTreesQuartets$(string(integer(time()/1000))).txt")
    writeSummary && descData(d,"summaryTreesQuartets.txt")
    return d
end

# rename the function readInputData to make it more user-friendly
"""
    readTrees2CF(treefile)
    readTrees2CF(vector of trees)

Read trees in parenthetical format from a file, or take a vector of trees already read,
and calculate the proportion of these trees having a given quartet (concordance factor: CF),
for all quartets or for a sample of quartets.
Optional arguments include:

- quartetfile: name of text file with list of 4-taxon subsets to be analyzed. If none is specified, the function will list all possible 4-taxon subsets.
- whichQ="rand": to choose a random sample of 4-taxon subsets
- numQ: size of random sample (ignored if whichQ is not set to "rand")
- writeTab=false: does not write the observedCF to a table (default true)
- CFfile: name of file to save the observedCF (default tableCF.txt)
- writeQ=true: save intermediate files with the list of all 4-taxon subsets and chosen random sample (default false).
- writeSummary: write descriptive stats of input data (default: true)
- nexus: if true, it assumes the gene trees are written in nexus file (default: false)
"""
function readTrees2CF(treefile::AbstractString; quartetfile="none"::AbstractString, whichQ="all"::AbstractString, numQ=0::Integer,
                      writeTab=true::Bool, CFfile="none"::AbstractString,
                      taxa=Vector{String}()::Union{Vector{String},Vector{Int}},
                      writeQ=false::Bool, writeSummary=true::Bool, nexus=false::Bool)
    if(nexus)
        trees = readNexusTrees(treefile)
    else
        trees = readInputTrees(treefile)
    end
    if length(taxa)==0        # unionTaxa(trees) NOT default argument:
      taxa = unionTaxa(trees) # otherwise: tree file is read twice
    end
    readTrees2CF(trees, quartetfile=quartetfile, whichQ=whichQ, numQ=numQ, writeTab=writeTab,
                 CFfile=CFfile, taxa=taxa, writeQ=writeQ, writeSummary=writeSummary)
end

# same as before, but with input vector of HybridNetworks
function readTrees2CF(trees::Vector{HybridNetwork}; quartetfile="none"::AbstractString, whichQ="all"::AbstractString, numQ=0::Integer, writeTab=true::Bool, CFfile="none"::AbstractString, taxa=unionTaxa(trees)::Union{Vector{String},Vector{Int}}, writeQ=false::Bool, writeSummary=true::Bool)
    whichQ == "all" || whichQ == "rand" ||
        error("whichQ should be all or rand, not $(whichQ)")
    if(quartetfile == "none")
        readInputData(trees, Symbol(whichQ), numQ, taxa, writeTab, CFfile, writeQ, writeSummary)
    else
        readInputData(trees, quartetfile, Symbol(whichQ), numQ, writeTab, CFfile, writeQ, writeSummary)
    end
end

# ---------------------- descriptive stat for input data ----------------------------------

# function to check how taxa is represented in the input trees
function taxaTreesQuartets(trees::Vector{HybridNetwork}, quartets::Vector{Quartet},s::IO)
    taxaT = unionTaxa(trees)
    taxaQ = unionTaxa(quartets)
    dif = symdiff(taxaT,taxaQ)
    isempty(dif) ? write(s,"\n same taxa in gene trees and quartets: $(taxaT)\n") :
                   write(s,"\n $(length(dif)) different taxa found in gene trees and quartets. \n Taxa $(intersect(taxaT,dif)) in trees, not in quartets; and taxa $(intersect(taxaQ,dif)) in quartets, not in trees\n")
    u = union(taxaT,taxaQ)
    for taxon in u
        numT = taxonTrees(taxon,trees)
        #numQ = taxonQuartets(taxon,quartets)
        write(s,"Taxon $(taxon) appears in $(numT) input trees ($(round(100*numT/length(trees),2)) %)\n")  #and $(numQ) quartets ($(round(100*numQ/length(quartets),2)) %)\n")
    end
end

taxaTreesQuartets(trees::Vector{HybridNetwork}, quartets::Vector{Quartet}) = taxaTreesQuartets(trees, quartets, STDOUT)

# function that counts the number of trees in which taxon appears
function taxonTrees(taxon::AbstractString, trees::Vector{HybridNetwork})
    suma = 0
    for t in trees
        suma += in(taxon,t.names) ? 1 : 0
    end
    return suma
end

# function that counts the number of quartets in which taxon appears
function taxonQuartets(taxon::AbstractString, quartets::Vector{Quartet})
    suma = 0
    for q in quartets
        suma += in(taxon,q.taxon) ? 1 : 0
    end
    return suma
end

# function to create descriptive stat from input data, will save in stream sout
# which can be a file or STDOUT
# default: send to STDOUT
# pc: only 4-taxon subsets with percentage of gene trees less than pc will be printed (default 70%)
function descData(d::DataCF, sout::IO, pc::Float64)
    0<=pc<=1 || error("percentage of missing genes should be between 0,1, not: $(pc)")
    if(!isempty(d.tree))
        print(sout,"data consists of $(d.numTrees) gene trees and $(d.numQuartets) 4-taxon subsets\n")
        taxaTreesQuartets(d.tree,d.quartet,sout)
        print(sout,"----------------------------\n\n")
        print(sout,"will print below only the 4-taxon subsets with data from <= $(round((pc)*100,2))% genes\n")
        for q in d.quartet
            percent  = q.ngenes == -1 ? 0.0 : round(q.ngenes/d.numTrees*100,2)
            if(percent < pc)
                print(sout,"4-taxon subset $(q.taxon) obsCF constructed with $(q.ngenes) gene trees ($(percent)%)\n")
            end
        end
        print(sout,"----------------------------\n\n")
    else
        if(!isempty(d.quartet))
            print(sout,"data consists of $(d.numQuartets) 4-taxon subsets")
            taxa=unionTaxa(d.quartet)
            print(sout,"\nTaxa: $(taxa)\n")
            print(sout,"Number of Taxa: $(length(taxa))\n")
            numQ = binomial(length(taxa),4);
            print(sout,"Maximum number of 4-taxon subsets: $(numQ). Thus, $(round(100*d.numQuartets/numQ,2)) percent of 4-taxon subsets sampled\n")
        end
    end
end

function descData(d::DataCF, filename::AbstractString,pc::Float64)
    println("descriptive stat of input data printed to file $(filename)")
    s = open(filename, "w")
    descData(d,s,pc)
    close(s)
end

descData(d::DataCF, sout::IO) = descData(d, sout,0.7)
descData(d::DataCF) = descData(d, STDOUT,0.7)
descData(d::DataCF,pc::Float64) = descData(d, STDOUT,pc)
descData(d::DataCF, filename::AbstractString) = descData(d, filename,0.7)

"""
`summarizeDataCF(d::DataCF)`

function to summarize the information contained in a DataCF object. It has the following optional arguments:
- filename: if provided, the summary will be saved in the filename, not to screen
- pc (number between (0,1)): threshold of percentage of missing genes to identify 4-taxon subsets with fewer genes than the threshold
"""
function summarizeDataCF(d::DataCF; filename="none"::AbstractString, pc=0.7::Float64)
    0<=pc<=1 || error("percentage of missing genes should be between 0,1, not: $(pc)")
    if(filename == "none")
        descData(d,STDOUT,pc)
    else
        descData(d,filename,pc)
    end
end

# ------------------ read starting tree from astral and put branch lengths ------------------------------

# function to read the starting topology (can be tree/network)
# if updateBL=true, updates the branch lengths with the obsCF in d
# by default, updateBL=true
function readStartTop(file::AbstractString,d::DataCF,updateBL::Bool)
    net = readTopologyUpdate(file)
    if(updateBL)
        updateBL!(net,d)
    end
    return net
end

"""
`readStartTop(treefile,d::DataCF)`

function to read a tree in parenthetical format from text file treefile and a DataCF object to update the branch lengths according to the average observed CF for any given edge.
"""
readStartTop(file::AbstractString,d::DataCF) = readStartTop(file,d,true)
#readStartTop(file::AbstractString) = readStartTop(file,DataCF(),false) #not sure why we need this one

# function to update starting branch lengths for starting tree read from ASTRAL
# BL are updated as -log(3/2(1-mean(obsCF)))
# based on Cecile's chr4-species-tree-units.r
# input: starting tree, data after read table of obsCF
function updateBL!(net::HybridNetwork,d::DataCF)
    if(isTree(net))
        parts = edgesParts(net)
        df = makeTable(net,parts,d)
        x=by(df,[:edge],df->DataFrame(Nquartets=length(df[:CF]),edgeL=-log(3/2*(1-mean(df[:CF])))))
        # ommitting columns: meanCF=mean(df[:CF]), sdCF=std(df[:CF])
        edges = x[:edge]
        lengths = x[:edgeL]
        for i in 1:length(edges)
            try
                ind = getIndexEdge(edges[i],net)
            catch
                error("edge $(edges[i]) not in net")
            end
            ind = getIndexEdge(edges[i],net)
            if (net.edge[ind].length < 0.0 || net.edge[ind].length==1.0)
                # readTopologyLevel1 changes missing branch length to 1.0
                setLength!(net.edge[ind], (lengths[i] > 0 ? lengths[i] : 0.0))
            end
        end
        for e in net.edge
            if e.length < 0.0 # some edges might have *no* quartet in the data
                setLength!(e, 1.0)
            end
        end
        return x
    else
        warn("updateStartBL was created for a tree, and net here is not a tree, so no branch lengths updated")
    end
end

# function to get part1,part2,part3,part4 for each edge in net.edge
# returns a EdgeParts object
function edgesParts(net::HybridNetwork)
    parts = EdgeParts[] #vector to hold part1,...,part4 for each edge
    for e in net.edge
        if(isInternalEdge(e))
            length(e.node) == 2 || error("strange edge with $(length(e.node)) nodes instead of 2")
            n1 = e.node[1]
            n2 = e.node[2]
            e11,e12 = hybridEdges(n1,e)
            e21,e22 = hybridEdges(n2,e)
            part1 = Node[]
            part2 = Node[]
            part3 = Node[]
            part4 = Node[]
            getDescendants!(getOtherNode(e11,n1),e11,part1)
            getDescendants!(getOtherNode(e12,n1),e12,part2)
            getDescendants!(getOtherNode(e21,n2),e21,part3)
            getDescendants!(getOtherNode(e22,n2),e22,part4)
            push!(parts, EdgeParts(e.number,part1,part2,part3,part4))
        end
    end
    return parts
end

# aux function to traverse the network from a node and an edge
# based on traverseContainRoot
# warning: it does not go accross hybrid node, minor hybrid edge
# there is another getDescendants in update.jl for updatePartition
function getDescendants!(node::Node, edge::Edge, descendants::Array{Node,1})
    if(node.leaf)
        push!(descendants, node)
    else
        for e in node.edge
            if(!isEqual(edge,e) && e.isMajor)
                other = getOtherNode(e,node);
                getDescendants!(other,e, descendants);
            end
        end
    end
end

# function to make table to later use in updateBL
# uses vector parts obtained from edgeParts function
function makeTable(net::HybridNetwork, parts::Vector{EdgeParts},d::DataCF)
    df = DataFrame(edge=Int[],t1=AbstractString[],t2=AbstractString[],t3=AbstractString[],t4=AbstractString[],resolution=AbstractString[],CF=Float64[])
    sortedDataQ = [sort(q.taxon) for q in d.quartet]
    for p in parts #go over internal edges too
        for t1 in p.part1
            for t2 in p.part2
                for t3 in p.part3
                    for t4 in p.part4
                        tx1 = net.names[t1.number]
                        tx2 = net.names[t2.number]
                        tx3 = net.names[t3.number]
                        tx4 = net.names[t4.number]
                        nam = [tx1,tx2,tx3,tx4]
                        snam = sort(nam)
                        # row = getIndex(true,[sort(nam) == sort(q.taxon) for q in d.quartet])
                        # getIndex was getting the first index only: like findfirst
                        row = findin([dnam==snam for dnam in sortedDataQ], true)
                        for r in row # nothing if tax set not found: length(row)=0
                          col,res = resolution(nam,d.quartet[r].taxon)
                          push!(df, [p.edgenum,tx1,tx2,tx3,tx4,res,d.quartet[r].obsCF[col]])
                        end
                    end
                end
            end
        end
    end
    return df
end

# function to determine the resolution of taxa picked from part1,2,3,4 and DataCF
# names: taxa from part1,2,3,4
# rownames: taxa from table of obsCF
function resolution(names::Vector{String},rownames::Vector{String})
    length(names) == length(rownames) || error("names and rownames should have the same length")
    length(names) == 4 || error("names should have 4 entries, not $(length(names))")
    bin = [n == names[1] || n == names[2] ? 1 : 0 for n in rownames]
    if(bin == [1,1,0,0] || bin == [0,0,1,1])
        return 1,"12|34"
    elseif(bin == [1,0,1,0] || bin == [0,1,0,1])
        return 2,"13|24"
    elseif(bin == [1,0,0,1] || bin == [0,1,1,0])
        return 3,"14|23"
    else
        error("strange resolution $(bin)")
    end
end

# function to extract a quartet from a matrix M
# obtained from tree2Matrix (defined in file compareNetworks.jl)
# this function is meant to replace extractQuartet! in calculateObsCFAll
# input: Quartet, Matrix, vector of taxa names
# returns 1 if quartet found is 12|34, 2 if 13|24, 3 if 14|23, and 0 if not found
function extractQuartetTree(q::Quartet, M::Matrix{Int},S::Union{Vector{String},Vector{Int}})
    DEBUG && println("extractQuartet: $(q.taxon)")
    DEBUG && println("matrix: $(M)")
    try
        ind1 = getIndex(q.taxon[1],S)
        ind2 = getIndex(q.taxon[2],S)
        ind3 = getIndex(q.taxon[3],S)
        ind4 = getIndex(q.taxon[4],S)
    catch
        error("some taxon in quartet $(q.taxon) not found in list of all species $(S)")
    end
    ind1 = getIndex(q.taxon[1],S)
    ind2 = getIndex(q.taxon[2],S)
    ind3 = getIndex(q.taxon[3],S)
    ind4 = getIndex(q.taxon[4],S)
    subM = M[:,[ind1+1,ind2+1,ind3+1,ind4+1]]
    DEBUG && println("subM: $(subM)")
    for r in 1:size(subM,1) #rows in subM
        DEBUG && println("subM[r,:]: $(subM[r,:])")
        if(subM[r,:] == [0,0,1,1] || subM[r,:] == [1,1,0,0])
            return 1
        elseif(subM[r,:] == [0,1,0,1] || subM[r,:] == [1,0,1,0])
            return 2
        elseif(subM[r,:] == [0,1,1,0] || subM[r,:] == [1,0,0,1])
            return 3
        end
    end
    return 0
end

# function that will give the qth quartet without making a list of all quartets
# input: n number of taxa, q desired index of quartet
# returns vector of int, e.g. 1234
function whichQuartet(n::Int, q::Int)
    p = 4
    q <= binom(n,p) || error("the index for the quartet $(q) needs to be less than choose(n,4)=$(binom(n,p))")
    n > 4 || error("there must be at least 5 taxa, not $(n)")
    quartet = Int[]
    while(n > 1)
        abs = binom(n-1,p) #fixit: we don't want to compute this, we want to look for it in a table
        if(q > abs)
            push!(quartet,n)
            n -= 1
            p -= 1
            q = q-abs
        else
            n -= 1
        end
    end
    if(length(quartet) == 3)
        push!(quartet,1)
    end
    quartet = quartet[[4,3,2,1]] #sort
    return quartet
end

# function to write a quartet on integer to taxon names
# it creates a Quartet type
# input: list of taxa names, vector of integers (each integer corresponds to a taxon name),
# num is the number of the Quartet
# assumes taxa is sorted already
function createQuartet(taxa::Union{Vector{String},Vector{Int}},qvec::Vector{Int}, num::Int)
    length(qvec) == 4 || error("a quartet should have only 4 taxa, not $(length(qvec))")
    names = String[]
    for i in qvec
        i <= length(taxa) || error("want taxon number $(i) in list of taxon names $(taxa) which has only $(length(taxa)) names")
        push!(names,string(taxa[i]))
    end
    return Quartet(num,names,[1.0,0.0,0.0])
end

## internal function to read a treefile in nexus format
function readNexusTrees(file::AbstractString)
    try
        s = open(file)
    catch
        error("Could not find or open $(file) file");
    end
    vnet = HybridNetwork[]
    s = open(file)
    numl = 1
    for line in eachline(s)
        line = strip(line) # remove spaces
        DEBUG && println("$(line)")
        c = isempty(line) ? "" : line[1:4]
        if(c == "Tree")
            init = search(line,'(')
            endit = search(line,';')
            line = line[init:endit]
           try
               push!(vnet, readTopologyUpdate(line,false))
               ##@show line
           catch(err)
               error("could not read tree in line $(numl). The error is $(err)")
           end
        end
        numl += 1
    end
    close(s)
    return vnet # consistent output type: HybridNetwork vector. might be of length 0.
end


text:
using Oceananigans: Flat, Bounded, RegularRectilinearGrid

ENGAGING_LESBRARY_DIR = "https://engaging-web.mit.edu/~alir/lesbrary/neural_free_convection_training_data"

SIMULATION_IDS = 1:9

LESBRARY_DATA_DEPS = [
    DataDep("free_convection_$id",
            "LESbrary.jl free convection simulation $id",
            joinpath(ENGAGING_LESBRARY_DIR, "free_convection_$id", "instantaneous_statistics_with_halos.jld2")
    ) for id in SIMULATION_IDS
]

function validate_simulation_ids(ids_train, ids_test)
    @info "Simulation IDs: $(collect(SIMULATION_IDS))"
    @info "Training simulations: $(collect(ids_train))"
    @info "Testing simulations: $(collect(ids_test))"

    ids_not_used = setdiff(SIMULATION_IDS, ids_train, ids_test)
    if !isempty(ids_not_used)
        @warn "Simulations not used: $(collect(ids_not_used))"
    end

    ids_intersection = intersect(ids_train, ids_test)
    if !isempty(ids_intersection)
        @warn "Simulations used for both training and testing: $(collect(ids_intersection))"
    end

    return nothing
end

function load_data(ids_train, ids_test, Nz)

    @info "Constructing FieldTimeSeries..."

    datasets = Dict{Int, FieldDataset}(
        id => FieldDataset(@datadep_str "free_convection_$id/instantaneous_statistics_with_halos.jld2"; ArrayType=Array{Float32}, metadata_paths=["parameters"])
        for id in SIMULATION_IDS
    )

    @info "Injecting surface fluxes..."

    for id in SIMULATION_IDS
        add_surface_fluxes!(datasets[id])
    end

    @info "Coarsening grid..."

    les_grid = datasets[1]["T"].grid

    topo = (Flat, Flat, Bounded)
    domain = (les_grid.zF[1], les_grid.zF[les_grid.Nz+1])
    coarse_grid = RegularRectilinearGrid(Float32, topology=topo, size=Nz, z=domain)

    @info "Coarse graining data..."

    coarse_datasets = Dict{Int, FieldDataset}(
        id => coarse_grain(ds, coarse_grid)
        for (id, ds) in datasets
    )

    @info "Partitioning data into training and testing datasets..."

    training_datasets = Dict{Int, FieldDataset}(id => datasets[id] for id in ids_train)
    testing_datasets = Dict{Int, FieldDataset}(id => datasets[id] for id in ids_test)

    coarse_training_datasets = Dict{Int, FieldDataset}(id => coarse_datasets[id] for id in ids_train)
    coarse_testing_datasets = Dict{Int, FieldDataset}(id => coarse_datasets[id] for id in ids_test)

    return (; datasets, coarse_datasets, training_datasets, testing_datasets, coarse_training_datasets, coarse_testing_datasets)
end


text:
struct ReverseCompoundPeriod <: AbstractTime
   cperiod::CompoundPeriod
end

reverse(x::CompoundPeriod) = ReverseCompoundPeriod(x)
reverse(x::ReverseCompoundPeriod) = x.cperiod
reverse(x::Period) = reverse(CompoundPeriod(x))

function string(rperiod::ReverseCompoundPeriod)
    join(([string(aperiod) for aperiod in rperiod]...,),", ")
end

function show(io::IO, x::ReverseCompoundPeriod)
   print(io, string(x))
end

const StdOut = isdefined(Base, :stdout) ? Base.stdout : Base.STDOUT

function show(x::ReverseCompoundPeriod)
    print(StdOut, string(x))
end


text:
using Formatting

function writeabq(filename::String, mesh::Mesh)
    # open file
    open(filename, "w") do afile
        # node header
        write(afile, "*Node\n")

        # node coordinates
        for (nodeid, nodecoords) in enumerate(mesh.nodes)
            write(afile, format("{}, {}, {}\n", nodeid, nodecoords...))
        end

        # elements
        for (elem_type, element_ids) in mesh.elementsets
            # element header
            write(afile, "*Element, type=$elem_type\n")
            for element_id in element_ids
                element = mesh.elements[element_id]
                write(afile, format("{}, {}, {}, {}, {}\n", element_id, element...))
            end
        end
    end
end


text:
function numback_hesvec!(dy, f, x, v, cache1 = similar(v), cache2 = similar(v))
    g = let f=f
        (dx, x) -> dx .= first(Zygote.gradient(f,x))
    end
    T = eltype(x)
    # Should it be min? max? mean?
    ϵ = sqrt(eps(real(T))) * max(one(real(T)), abs(norm(x)))
    @. x += ϵ*v
    g(cache1,x)
    @. x -= 2ϵ*v
    g(cache2,x)
    @. dy = (cache1 - cache2)/(2ϵ)
end

function numback_hesvec(f, x, v)
    g = x -> first(Zygote.gradient(f,x))
    T = eltype(x)
    # Should it be min? max? mean?
    ϵ = sqrt(eps(real(T))) * max(one(real(T)), abs(norm(x)))
    x += ϵ*v
    gxp = g(x)
    x -= 2ϵ*v
    gxm = g(x)
    (gxp - gxm)/(2ϵ)
end

function autoback_hesvec!(dy, f, x, v,
                          cache2 = Dual{typeof(ForwardDiff.Tag(DeivVecTag,eltype(x))),eltype(x),1}.(x, ForwardDiff.Partials.(Tuple.(reshape(v, size(x))))),
                          cache3 = Dual{typeof(ForwardDiff.Tag(DeivVecTag,eltype(x))),eltype(x),1}.(x, ForwardDiff.Partials.(Tuple.(reshape(v, size(x))))))
    g = let f=f
        (dx, x) -> dx .= first(Zygote.gradient(f,x))
    end
    cache2 .= Dual{typeof(ForwardDiff.Tag(DeivVecTag,eltype(x))),eltype(x),1}.(x, ForwardDiff.Partials.(Tuple.(reshape(v, size(x)))))
    g(cache3,cache2)
    dy .= partials.(cache3, 1)
end

function autoback_hesvec(f, x, v)
    g = x -> first(Zygote.gradient(f,x))
    y = Dual{typeof(ForwardDiff.Tag(DeivVecTag,eltype(x))),eltype(x),1}.(x, ForwardDiff.Partials.(Tuple.(reshape(v, size(x)))))
    ForwardDiff.partials.(g(y), 1)
end


text:
struct NoFeasibleInfill <: Exception
    msg::String
end

mutable struct Acquire
    global_opt::KrigingModel.Optimizer
    local_opt::KrigingModel.Optimizer
end

function Acquire(;
    global_optimizer=ISRESOptimizer(; maxeval=10000, ftol_rel=1e-2),
    local_optimizer=LBFGSOptimizer(; maxeval=1000),
)
    return Acquire(global_optimizer, local_optimizer)
end

function doacquire(opt::Acquire, obj, lb::AbstractVector, ub::AbstractVector)
    fx, x, ret = maximize(opt.global_opt, obj, :Auto, lb, ub, (lb + ub) / 2)
    fx, x, ret = maximize(opt.local_opt, obj, :Auto, lb, ub, x)
    return x, fx
end

function acquire(
    krg::Kriging,
    a::EGOAquisition,
    opt::Acquire,
    lb::AbstractVector,
    ub::AbstractVector,
    x_cache::AbstractVector...,
)
    f = x -> begin
        y = objective(a, krg, x)
        for xc in x_cache
            y *= 1.0 - exp(-sum((x - xc) .^ 2))
        end
        y
    end
    x, fx = doacquire(opt, f, lb, ub)
    if fx == 0
        throw(NoFeasibleInfill("Cannot find feasible infill. Kriging model is converged."))
    end
    return x, fx
end


text:
module NCTiles

using NCDatasets, NetCDF, MeshArrays
using Dates, Printf, Pkg, LazyArtifacts
#using MITgcmTools

version()=Pkg.TOML.parsefile(joinpath(dirname(pathof(NCTiles)), "..", "Project.toml"))["version"]

p=dirname(pathof(NCTiles))
artifact_toml = joinpath(p, "../Artifacts.toml")
NCTILES_TESTCASES_hash = artifact_hash("NCTILES_TESTCASES", artifact_toml)
NCTILES_TESTCASES = joinpath(artifact_path(NCTILES_TESTCASES_hash)*"/","nctiles-testcases-0.1/")

"""
    pkg_version(name::String)

Get the version number (in current environment) for a package name
"""
function pkg_version(name::String)
    tmp1=Pkg.dependencies()
    tmp2=Vector{Union{Missing, VersionNumber}}(missing, 1)
    for i in collect(keys(tmp1))
        if tmp1[i].name==name
            tmp2[1]=tmp1[i].version
        end
    end
    return tmp2[1]
end

include("write.jl")
include("read.jl")
include("tile_support.jl")
include("helper_functions.jl")

export NCT, NCvar, BinData, NCData, TileData, readbin, readncfile
export createfile, addDim, addVar, addData, addDimData
export write

end # module


text:
@testset "Solutions" begin
  @testset "Estimation" begin
    d = PointSet(rand(2,100))
    m = Dict(:z => 1:100)
    v = Dict(:z => 1:100)
    s = EstimationSolution(d, m, v)
    @test s[:z] == (mean=1:100, variance=1:100)

    d = RegularGrid{Float64}(10,10)
    s = EstimationSolution(d, m, v)
    @test s[:z] == (mean=reshape(1:100,10,10), variance=reshape(1:100,10,10))

    @test sprint(show, s) == "2D EstimationSolution"
    @test sprint(show, MIME"text/plain"(), s) == "2D EstimationSolution\n  domain: 10×10 RegularGrid{Float64,2}\n  variables: z"

    if visualtests
      gr(size=(800,400))
      @plottest plot(s) joinpath(datadir,"EstimationSolution.png") !istravis
    end
  end

  @testset "Simulation" begin
    # TODO
  end

  @testset "Learning" begin
    # TODO
  end
end


text:
######################################################
# Created by Alex Dowling (alexdowling.net)
# while at the University of Wisconsin-Madison
# ####################################################

# Define problem
include("../DegeneracyHunter.jl")
using JuMP

# Define problem
function lp1()

	m = Model()

	@variable(m, 0 <= x[1:3] <= 5)
	@variable(m, 0 <= y <= 0)
	@constraint(m, x[1] + x[2] >= 1)
	@constraint(m, x[1] + x[2] + x[3] == 1)
	@constraint(m, x[2] - 2*x[3] <= 1)
	@constraint(m, x[1] + x[3] >= 1)
	@constraint(m, x[1] + x[2] + x[3] == 1) # Redundant constraint - makes problem degenerate
	@objective(m, Min, sum{x[i],i=1:3})

	return m
end

function initialize!(m::Model)
	x = getvariable(m,:x)
	setvalue(x[1], 1.0)
	setvalue(x[2], 5.0)
	setvalue(x[3], -1.0)

	# Intentionally do not initialize y

	return nothing
end

# Create and initialize model
m = lp1()
initialize!(m)

# Check initialization
DegeneracyHunter.printVariableDiagnostics(m)

# Print equations with absolute residuals greater than 0.001
DegeneracyHunter.printInfeasibleEquations(m, 0.001)

# Check for degeneracy constraints
DegeneracyHunter.degeneracyHunter(m, includeBounds=true)

# Solve model
tic()
status = solve(m)
tm = toq()

# Print inactive constraints at solution
DegeneracyHunter.printInactiveEquations(m)

# Package and print problem size information
ps = DegeneracyHunter.assembleProblemStats(m,status,tm)
DegeneracyHunter.printProblemStats(ps)


text:
function DSSTBatchLSModel(arg0::Vector{IntegratedPropagatorBuilder}, arg1::List, arg2::ParameterDriversList, arg3::ModelObserver, arg4::PropagationType, arg5::PropagationType)
    return DSSTBatchLSModel((Vector{IntegratedPropagatorBuilder}, List, ParameterDriversList, ModelObserver, PropagationType, PropagationType), arg0, arg1, arg2, arg3, arg4, arg5)
end

function create_propagators(obj::DSSTBatchLSModel, arg0::RealVector)
    return jcall(obj, "createPropagators", Vector{DSSTPropagator}, (RealVector,), arg0)
end

function fetch_evaluated_measurement(obj::DSSTBatchLSModel, arg0::jint, arg1::EstimatedMeasurement)
    return jcall(obj, "fetchEvaluatedMeasurement", void, (jint, EstimatedMeasurement), arg0, arg1)
end

function get_evaluations_count(obj::DSSTBatchLSModel)
    return jcall(obj, "getEvaluationsCount", jint, ())
end

function get_iterations_count(obj::DSSTBatchLSModel)
    return jcall(obj, "getIterationsCount", jint, ())
end

function get_selected_propagation_drivers_for_builder(obj::DSSTBatchLSModel, arg0::jint)
    return jcall(obj, "getSelectedPropagationDriversForBuilder", ParameterDriversList, (jint,), arg0)
end

function is_forward_propagation(obj::DSSTBatchLSModel)
    return jcall(obj, "isForwardPropagation", jboolean, ())
end

function set_evaluations_counter(obj::DSSTBatchLSModel, arg0::Incrementor)
    return jcall(obj, "setEvaluationsCounter", void, (Incrementor,), arg0)
end

function set_iterations_counter(obj::DSSTBatchLSModel, arg0::Incrementor)
    return jcall(obj, "setIterationsCounter", void, (Incrementor,), arg0)
end

function value(obj::DSSTBatchLSModel, arg0::RealVector)
    return jcall(obj, "value", Pair, (RealVector,), arg0)
end


text:
function with_temp_dir(f::Function)
    original_working_directory = pwd()

    temp_dir = mktempdir()
    atexit(() -> rm(temp_dir; force=true, recursive=true))

    cd(temp_dir)
    result = f(temp_dir)

    cd(original_working_directory)
    rm(temp_dir; force=true, recursive=true)
    return result
end

function with_temp_depot(f::Function)
    original_depot_path = deepcopy(Base.DEPOT_PATH)
    result = with_temp_dir() do temp_depot
        empty!(Base.DEPOT_PATH)
        push!(Base.DEPOT_PATH, temp_depot)
        return f()
    end
    empty!(Base.DEPOT_PATH)
    for x in original_depot_path
        push!(Base.DEPOT_PATH, x)
    end
    return result
end