text:
# This file is a part of Julia. License is MIT: https://julialang.org/license

"""
    ConsoleLogger(stream=stderr, min_level=Info; meta_formatter=default_metafmt,
                  show_limited=true, right_justify=0)

Logger with formatting optimized for readability in a text console, for example
interactive work with the Julia REPL.

Log levels less than `min_level` are filtered out.

Message formatting can be controlled by setting keyword arguments:

* `meta_formatter` is a function which takes the log event metadata
  `(level, _module, group, id, file, line)` and returns a color (as would be
  passed to printstyled), prefix and suffix for the log message.  The
  default is to prefix with the log level and a suffix containing the module,
  file and line location.
* `show_limited` limits the printing of large data structures to something
  which can fit on the screen by setting the `:limit` `IOContext` key during
  formatting.
* `right_justify` is the integer column which log metadata is right justified
  at. The default is zero (metadata goes on its own line).
"""
struct ConsoleLogger <: AbstractLogger
    stream::IO
    min_level::LogLevel
    meta_formatter
    show_limited::Bool
    right_justify::Int
    message_limits::Dict{Any,Int}
end
function ConsoleLogger(stream::IO=stderr, min_level=Info;
                       meta_formatter=default_metafmt, show_limited=true,
                       right_justify=0)
    ConsoleLogger(stream, min_level, meta_formatter,
                  show_limited, right_justify, Dict{Any,Int}())
end

shouldlog(logger::ConsoleLogger, level, _module, group, id) =
    get(logger.message_limits, id, 1) > 0

min_enabled_level(logger::ConsoleLogger) = logger.min_level

# Formatting of values in key value pairs
showvalue(io, msg) = show(io, "text/plain", msg)
function showvalue(io, e::Tuple{Exception,Any})
    ex,bt = e
    showerror(io, ex, bt; backtrace = bt!=nothing)
end
showvalue(io, ex::Exception) = showerror(io, ex)

function default_logcolor(level)
    level < Info  ? Base.debug_color() :
    level < Warn  ? Base.info_color()  :
    level < Error ? Base.warn_color()  :
                    Base.error_color()
end

function default_metafmt(level, _module, group, id, file, line)
    color = default_logcolor(level)
    prefix = (level == Warn ? "Warning" : string(level))*':'
    mod = _module === nothing ? "" : "$(_module) "
    suffix = (Info <= level < Warn) ? "" : "@ $mod$(basename(file)):$line"
    color,prefix,suffix
end

# Length of a string as it will appear in the terminal (after ANSI color codes
# are removed)
function termlength(str)
    N = 0
    in_esc = false
    for c in str
        if in_esc
            if c == 'm'
                in_esc = false
            end
        else
            if c == '\e'
                in_esc = true
            else
                N += 1
            end
        end
    end
    return N
end

function handle_message(logger::ConsoleLogger, level, message, _module, group, id,
                        filepath, line; maxlog=nothing, kwargs...)
    if maxlog != nothing && maxlog isa Integer
        remaining = get!(logger.message_limits, id, maxlog)
        logger.message_limits[id] = remaining - 1
        remaining > 0 || return
    end

    # Generate a text representation of the message and all key value pairs,
    # split into lines.
    msglines = [(indent=0,msg=l) for l in split(chomp(string(message)), '\n')]
    dsize = displaysize(logger.stream)
    if !isempty(kwargs)
        valbuf = IOBuffer()
        rows_per_value = max(1, dsize[1]÷(length(kwargs)+1))
        valio = IOContext(IOContext(valbuf, logger.stream),
                          :displaysize=>(rows_per_value,dsize[2]-5))
        if logger.show_limited
            valio = IOContext(valio, :limit=>true)
        end
        for (key,val) in pairs(kwargs)
            showvalue(valio, val)
            vallines = split(String(take!(valbuf)), '\n')
            if length(vallines) == 1
                push!(msglines, (indent=2,msg=SubString("$key = $(vallines[1])")))
            else
                push!(msglines, (indent=2,msg=SubString("$key =")))
                append!(msglines, ((indent=3,msg=line) for line in vallines))
            end
        end
    end

    # Format lines as text with appropriate indentation and with a box
    # decoration on the left.
    color,prefix,suffix = logger.meta_formatter(level, _module, group, id, filepath, line)
    minsuffixpad = 2
    buf = IOBuffer()
    iob = IOContext(buf, logger.stream)
    nonpadwidth = 2 + (isempty(prefix) || length(msglines) > 1 ? 0 : length(prefix)+1) +
                  msglines[end].indent + termlength(msglines[end].msg) +
                  (isempty(suffix) ? 0 : length(suffix)+minsuffixpad)
    justify_width = min(logger.right_justify, dsize[2])
    if nonpadwidth > justify_width && !isempty(suffix)
        push!(msglines, (indent=0,msg=SubString("")))
        minsuffixpad = 0
        nonpadwidth = 2 + length(suffix)
    end
    for (i,(indent,msg)) in enumerate(msglines)
        boxstr = length(msglines) == 1 ? "[ " :
                 i == 1                ? "┌ " :
                 i < length(msglines)  ? "│ " :
                                         "└ "
        printstyled(iob, boxstr, bold=true, color=color)
        if i == 1 && !isempty(prefix)
            printstyled(iob, prefix, " ", bold=true, color=color)
        end
        print(iob, " "^indent, msg)
        if i == length(msglines) && !isempty(suffix)
            npad = max(0, justify_width - nonpadwidth) + minsuffixpad
            print(iob, " "^npad)
            printstyled(iob, suffix, color=:light_black)
        end
        println(iob)
    end

    write(logger.stream, take!(buf))
    nothing
end

max_stars_repo_name: helgee/julia

max_issues_repo_name: helgee/julia

max_forks_repo_name: helgee/julia


text:
macro _primitive_columns(args...)
    funcs = Expr[]
    for arg in args
        arg_string = string(arg)
        push!(funcs, quote is_ch_type(::Val{Symbol($arg_string)}) = true end)
        push!(funcs, quote result_type(::Val{Symbol($arg_string)}) =
                                                            Vector{$arg} end)
        push!(funcs, quote
            function read_col_data(sock::ClickHouseSock,
                     num_rows::VarUInt,
                    ::Val{Symbol($arg_string)})
                return chread(sock, Vector{$arg}, num_rows)
            end
        end )
        push!(funcs, quote
            function write_col_data(sock::ClickHouseSock,
                data::AbstractVector{$arg},
                ::Val{Symbol($arg_string)})
                return chwrite(sock, data)
            end
        end )
        push!(funcs, quote
            function write_col_data(sock::ClickHouseSock,
                data::AbstractVector,
                ::Val{Symbol($arg_string)})
                return chwrite(sock, convert(Vector{$arg},data))
            end
        end )
        push!(funcs, quote deserialize(::Val{Symbol($arg_string)}) = $arg end )
    end
    return esc(:($(funcs...),))
end

@_primitive_columns(
    UInt8,
    UInt16,
    UInt32,
    UInt64,
    Int8,
    Int16,
    Int32,
    Int64,
    Float32,
    Float64,
    String
)

ext: jl


text:
"""
At the beginning of a page that uses a Unicode character encoding you may find some bytes
that represent the Unicode code point U+FEFF BYTE ORDER MARK (abbreviated as BOM ).

This helper function just skips BOM if present.
"""
function consumeBOM!(io)
    # BOM character detection
    startpos = position(io)
    if !eof(io) && Parsers.peekbyte(io) == 0xef
        Parsers.readbyte(io)
        (!eof(io) && Parsers.readbyte(io) == 0xbb) || Parsers.fastseek!(io, startpos)
        (!eof(io) && Parsers.readbyte(io) == 0xbf) || Parsers.fastseek!(io, startpos)
    end
    return
end

const ASCII_SPACE = UInt8(' ')
const ASCII_TAB = UInt8('\t')
const ASCII_RETURN = UInt8('\r')
const ASCII_NEWLINE = UInt8('\n')

"""
    ws = *(
            %x20 /              ; Space
            %x09 /              ; Horizontal tab
            %x0A /              ; Line feed or New line
            %x0D )              ; Carriage return
"""
@inline iswhitespace(c) = c == ASCII_SPACE  ||
                          c == ASCII_TAB ||
                          c == ASCII_RETURN ||
                          c == ASCII_NEWLINE
"""
`resize!` will be called rarely and it's better not to be inlined to avoid code explosion
"""
@noinline resize_vec!(vec::AbstractVector{T}, new_size) where {T} = resize!(vec, new_size)


text:
module MagSim

using CSV
using DataFrames
using Random
using Statistics

import Base: show

export

    # AbstractModels
    MagnetismModels,
    Ising,
    Potts,
    XY,

    # analysis_methods
    sim_dict_to_df,

    # cluster_algorithms
    wolff!,
    swendsen_wang!,

    # local_algorithms
    metropolis!,
    heat_bath!,

    # parallel_simulation_methods
    simulate_Ising_parallel,
    simulate_Potts_parallel,

    # plotting_methods
    plot_energy,
    plot_heat_capacity,
    plot_magnetization,
    plot_susceptibility,
    plot_Ising_lattice

# includes
include("./Observables.jl")
include("./Parameters.jl")
include("./MagnetismModels.jl")
include("./analysis_methods.jl")
include("./cluster_algorithms.jl")
include("./initialization_methods.jl")
include("./local_algorithms.jl")
include("./misc_methods.jl")
include("./parallel_simulation_methods.jl")
include("./plotting_methods.jl")
include("./update_methods.jl")

end # module


text:
function playEpisode(env::AbsEnvironment, agent::AbsAgent; verbose=false, randomInitial = false, learn=true, threshold=100000, render=false)
	state = randomInitial ? rand(getAllStates(env)) : getInitialState(env)
	totalRewards = 0.0
	numOfStates = 1

	verbose && println("Initial State: $(state)")

	while !isTerminal(state, env) && numOfStates < threshold
		render && render_env(env)
		action = play(agent, state, env; learn=learn)
		verbose && println("Action: $(action)")
		state, reward = transfer(env, state, action)
		verbose && println("State: $(state)")
		verbose && println("Reward: $(reward)")
		observe(agent, state, reward, env; learn=learn, terminal = (numOfStates + 1 == threshold || isTerminal(state, env)))
		totalRewards += reward
		numOfStates += 1
	end

	return (totalRewards, numOfStates)
end

max_stars_repo_name: ozanarkancan/ReinforcementLearning

max_issues_repo_name: ozanarkancan/ReinforcementLearning

max_forks_repo_name: ozanarkancan/ReinforcementLearning


text:
##########################
# Approximation problems
##########################
module ApproximationProblems
using ..Platforms
import ..Platforms: SamplingStyle, SolverStyle, dictionary, platform
import BasisFunctions: coefficienttype, element, elements
using BasisFunctions: Dictionary, TensorProductDict
using ..ExtensionFrames
using DomainSets: Domain

export ApproximationProblem
"""
    abstract type ApproximationProblem end

Approximation problem types group the arguments that a user supplies
to the `Fun` constructor. This may be a dictionary or a platform, for example.

A concrete approximaton problem implements an interface. For example, one can
ask for its sampling operator. However, all functionality is delegated to the
underlying dictionary or platform. Thus, an approximation problem is nothing but
an empty intermediate layer that passes on the questions to the right object.

In turn, by default a platform delegates all queries to the dictionary it
represents. Hence, all defaults in the interface are really specified by the
dictionary. The platform can intercept and change anything, for example in order
to specify a sampling operator that is different from the default of the dictionary.

Other routines besides `Fun` can group their arguments into an approximation
problem too. This ensures that these routines automatically implement the exact
same interface as the `Fun` constructor. Thus, it is easy to write a
`samplingoperator(...)` routine that accepts the exact same arguments as `Fun`,
and that simply returns the sampling operator that the `Fun` constructor would
use to solve the approximation problem.
"""
abstract type ApproximationProblem end

export samplingparam
"""
    samplingparam(ap::ApproximationProblem)

Return the sampling parameter if it was saved by the approximation problem.
Otherwise return `nothing`
"""
samplingparam(ap::ApproximationProblem) = ap.samplingparam

export setsamplingparam!
"""
    setsamplingparam!(ap::ApproximationProblem, L)

Save the sampling parameter in the approximation problem
"""
setsamplingparam!(ap::ApproximationProblem, L) = (ap.samplingparam = L)

export AbstractPlatformApproximation

"""
    abstract type AbstractPlatformApproximation <: ApproximationProblem

An approximation problem with a `platform` and a `parameter`
"""
abstract type AbstractPlatformApproximation <: ApproximationProblem end
platform(ap::AbstractPlatformApproximation) = ap.platform
export parameter
parameter(ap::AbstractPlatformApproximation) = ap.param
coefficienttype(ap::AbstractPlatformApproximation) = coefficienttype(ap.dict)

export PlatformApproximation
"""
    mutable struct PlatformApproximation <: ApproximationProblem

A `PlatformApproximation` stores a platform and a parameter value. This
corresponds to a specific dictionary (which is also stored). All queries are
delegated to the platform, and thus the answers may differ from the defaults.

As with dictionary approximations, the settings can be overruled by explicitly
specifying optional arguments to `Fun`.
"""
mutable struct PlatformApproximation <: AbstractPlatformApproximation
    platform    ::  Platform
    dict        ::  Dictionary
    param
    samplingparam

    PlatformApproximation(platform, param, samplingparam = nothing) =
        new(platform, dictionary(platform, param), param, samplingparam)
end

dictionary(ap::PlatformApproximation) = ap.dict

SamplingStyle(ap::PlatformApproximation) = SamplingStyle(ap.platform)
SolverStyle(samplingstyle::SamplingStyle, ap::PlatformApproximation) = SolverStyle(ap.platform, samplingstyle)

export AdaptiveApproximation
"""
    struct AdaptiveApproximation <: ApproximationProblem

An `AdaptiveApproximation` only stores a platform, with which to compute adaptive
approximations.
"""
struct AdaptiveApproximation <: ApproximationProblem
    platform    ::  Platform
end

dictionary(ap::AdaptiveApproximation, param) = ap.platform[param]

export DictionaryApproximation
mutable struct DictionaryApproximation <: ApproximationProblem
end

export ProductPlatformApproximation
"""
    mutable struct ProductPlatformApproximation{N} <: ApproximationProblem

A `ProductPlatformApproximation` corresponds to a product platform.
"""
mutable struct ProductPlatformApproximation{N} <: AbstractPlatformApproximation
    platform        ::  ProductPlatform{N}
    param
    productparam    ::  NTuple{N,Any}
    dict            ::  TensorProductDict

    samplingparam

    ProductPlatformApproximation{N}(platform::ProductPlatform, param, samplingparam=nothing) where {N} =
        new(platform, param, productparameter(platform, param), dictionary(platform, param), samplingparam)
end

setsamplingparam!(ap::ProductPlatformApproximation, param) =
    (@assert length(ap.productparam)==length(param);ap.samplingparam = param)

approximationproblem(platform::ProductPlatform{N}, param) where {N} =
    ProductPlatformApproximation{N}(platform, param)

SamplingStyle(ap::ProductPlatformApproximation) = SamplingStyle(ap.platform)
SolverStyle(samplingstyle::SamplingStyle, ap::ProductPlatformApproximation) = SolverStyle(ap.platform, samplingstyle)

elements(ap::ProductPlatformApproximation) = nothing == samplingparam(ap) ?
    error("Not possible to get elements of `ProductPlatformApproximation` without `samplingparam`") :
    length(samplingparam(ap))==length(ap.productparam) ?
    map(approximationproblem, elements(ap.platform), ap.productparam, samplingparam(ap)) :
    error("sampling parameter should contain $(length(ap.productparam)) elements")

unsafe_elements(ap::ProductPlatformApproximation) = map(approximationproblem, elements(ap.platform), ap.productparam)

dictionary(ap::ProductPlatformApproximation) = ap.dict

export approximationproblem
"""
    approximationproblem(args...)

# examples

## Create an approxmationproblem from a dictionary
```jldocs
julia> approximationproblem(Fourier(10))
```

## Create an approxmationproblem from a dictionary and a domain (extensionframe)
```jldocs
julia> approximationproblem(Fourier(10),0.0..0.5)
```

## Create an approxmationproblem from a platform
```jldocs
julia> approximationproblem(FourierPlatform())
```

## Create an approxmationproblem from a platform and a platform parameter
```jldocs
julia> approximationproblem(FourierPlatform(), 10)
```
"""
approximationproblem(dict::Dictionary) =
    approximationproblem(platform(dict), param(dict))
# This two-argument routine allows to specify a coefficient type, in which case the
# dictionary will be promoted if necessary. This allows the `Fun` constructor to ensure
# that the dictionary can handle complex coefficients, for example.
# approximationproblem(::Type{T}, dict::Dictionary) where {T} =
#     DictionaryApproximation(ensure_coefficienttype(T, dict))

# 2. An approximation problem from a platform and a concrete parameter value
approximationproblem(platform::Platform, param) = (@assert correctparamformat(platform, param);PlatformApproximation(platform, param))
approximationproblem(platform::Platform, param, L) = (@assert correctparamformat(platform, param);PlatformApproximation(platform, param, L))

# 3. An adaptive approximation problem from a platform
approximationproblem(platform::Platform) = AdaptiveApproximation(platform)

# From the adaptive approximation problem we can create a concrete platform
# approximation by specifying a parameter value.
approximationproblem(ap::AdaptiveApproximation, param) = approximationproblem(ap.platform, param)
approximationproblem(ap::AdaptiveApproximation, param, L) = approximationproblem(ap.platform, param, L)
end


text:
using Utilities,StaticNets,JLD
using PyCall; pygui(:qt); using PyPlot

## Load data
#whole period data file
loadFileName = "/home/Domenico/Dropbox/Dynamic_Networks/data/emid_data/juliaFiles/Weekly_eMid_Data_from_2009_06_22_to_2015_02_27.jld"
@load(loadFileName, AeMidWeekly_Active_T,banksIDs,inactiveBanks, YeMidWeekly_Active_T,weekInd,datesONeMid ,degsIO_T,strIO_T)

## Methods definitions

struct SnapSeqNetDirBin1
    obsT:: Array{<:Real,2}  # In degrees and out Degrees for each t . The last
                            # dimension is always time. 2NxT
    Par::Array{<:Real,2} # In and Out Parameters [InW_T;outW_T] for each t
end
SnapSeqNetDirBin1(obs::Array{<:Real,2}) = SnapSeqNetDirBin1(obs,zeros(obs))

function estimate(Model::SnapSeqNetDirBin1; degsIO_T::Array{<:Real,2}=Model.obsT,targetErr::Real=1e-12)
   N2,T = size(degsIO_T);N = round(Int,N2/2)
   UnParT = zeros(2N,T)
   prog = 0
   for t = 1:T
       degs_t = degsIO_T[:,t]
       UnParT[:,t]  = identify!(fooErgmDirBin1,estimate(fooErgmDirBin1; degIO = degs_t , targetErr =  targetErr)[1])
       round(t/T,2)>prog ? (prog=round(t/T,2);println(prog) ):()
   end
 return UnParT
end

function bigNegParFun(Model::ErgmDirBin1,N::Int; smallProb = 1e-1)
    #define a number big enough to play the role of Inf for
    # purposes of sampling N(N-1) bernoulli rvs with prob 1/(1+exp(  bigNumb))
    # prob of sampling degree > 0 (N-1)/(1+exp(bigNumb)) < 1e-6
    bigNumb = - log((N-1)/smallProb - 1)
    return bigNumb
 end

## Application to Data
snapModData = SnapSeqNetDirBin1(degsIO_T)
singSnapEst_T = estimate(snapModData)
N2,T = size(degsIO_T);N = round(Int,N2/2)

## Stima local likelihood, la ottengo facendo i locally weighted degrees
Nalpha = 30 #Number of different alpha values for EWMA
degsIO_EWMA_T_Alpha = zeros(Float64,2N,T,Nalpha)
alphaVals = linspace(0.4+1/Nalpha,1,Nalpha)
# compute the EWMA degrees for different values of alpha
for alphaInd=1:Nalpha
    alpha = alphaVals[alphaInd]
    alphaVec =ones(1) *alpha
    for t=1:T
        wVec = alphaVec./sum(alphaVec)
        alphaVec*=alpha; push!(alphaVec,alpha)
        #println(wVec)
        t>1? degsIO_EWMA_T_Alpha[:,t,alphaInd] = sumSq(wVec'.*degsIO_T[:,1:t],2):()# (degIO_EWMA_T[:,t] = degsIO_T[1,:])
    end
end
# estimate the single snapshots on the EWMA degrees
singSnapEst_T_Alpha = zeros(2N,T,Nalpha)
for n=1:Nalpha
    singSnapEst_T_Alpha[:,:,n] = estimate( SnapSeqNetDirBin1(degsIO_EWMA_T_Alpha[:,:,n]))
    println(n)
end

plot(singSnapEst_T[1:N,:]')
#starting in a uniform point do I always end up is some particular identification?
x = (sumSq(singSnapEst_T[1:N,:],1)-sumSq(singSnapEst_T[1+N:2N,:],1))
y = Float64.(sumSq(degsIO_T[:,:],1))
plot(x,y)
cor(x,y)

##
degIO = degsIO_T[:,1]
parUn,it = estimate(fooErgmDirBin1; degIO = degIO , targetErr =  1e-15)

#

max_stars_repo_name: domenicodigangi/SDERGM

max_issues_repo_name: domenicodigangi/SDERGM

max_forks_repo_name: domenicodigangi/SDERGM


text:
using Test, RDatasets
using JDF, Tables

iris = dataset("datasets", "iris")
ok = savejdf(iris, "plsdel.jdf")

Tables.columns(ok)
Tables.schema(ok)
propertynames(ok)
getproperty(ok, :Species)

rm("plsdel.jdf", force=true, recursive=true)


text:
#sin bars

using CairoMakie

f = Figure()
Axis(f[1, 1])

xs = 1:0.2:10
ys = 0.5 .* sin.(xs)

barplot!(xs, ys, color = :seagreen, strokecolor = :black, strokewidth = 1)
barplot!(xs, ys .- 1, fillto = -1, color = xs, strokecolor = :black, strokewidth = 1)

f

max_stars_repo_name: Ellon-M/visualizations

max_issues_repo_name: Ellon-M/visualizations

max_forks_repo_name: Ellon-M/visualizations


text:
# Parses distribution expression and returns expression for a segment cost function
# Required for PELT and BS macros
function cost_function(data::Any, dist_expr::Expr)
    if !Meta.isexpr(dist_expr,:call)
        error("Syntax error: expected distribution as argument")
    end

    dist_type = dist_expr.args[1]

    if dist_type == :Normal
        if length(dist_expr.args) != 3
            error("Normal distribution has two parameters")
        end
        μ, σ = dist_expr.args[2:3]
        if μ == :(:?) && σ != :(:?)
            println("Changepoint distribution is Normal with changing mean and fixed variance")
            return :(NormalMeanSegment($data , $σ))
        elseif μ != :(:?) && σ == :(:?)
            println("Changepoint distribution is Normal with fixed mean and changing variance")
            #μ = eval(Main, μ)
            return :(NormalVarSegment($data, $μ))
        elseif μ == :(:?) && σ == :(:?)
            println("Changepoint distribution is Normal with changing mean and changing variance")
            return :(NormalMeanVarSegment($data))
        else
            error("Must mark at least one Normal parameter as changing with a :? symbol")
        end

    elseif dist_type == :Exponential
        println("Changepoint distribution is Exponential with changing mean")
        return :(ExponentialSegment($data))

    elseif dist_type == :Poisson
        println("Changepoint distribution is Poisson with changing mean")
        return :(PoissonSegment($data))

    elseif dist_type == :Gamma
        if length(dist_expr.args) != 3
            error("Gamma distribution has two parameters")
        end
        alpha , beta = dist_expr.args[2:3]
        if alpha == :(:?) && beta != :(:?)
            println("Changepoint distribution is Gamma with changing shape and fixed rate")
            return :(GammaShapeSegment($data, $beta))
        elseif alpha != :(:?) && beta == :(:?)
            println("Changepoint distribution is Gamma with fixed shape and changing rate")
            return :(GammaRateSegment($data, $alpha))
        else
            error("Must mark at least one Gamma parameter as changing with a :? symbol")
        end

    elseif dist_type == :Nonparametric
        K = dist_expr.args[2]
        println("Changepoint method is Nonparametric")
        return :(NonparametricSegment($data, $K))

    elseif dist_type == :OLS
        println("Changepoints in peicewise linear regressions")
        return :(OLSSegment($data))

    else
        error("Distribution $(dist_type) has no implemented cost functions")
    end
end

"""
Creates a segment cost function given data and changepoint model expression

# Usage
@segement_cost data changepoint_model

# Cost functions

A changepoint model is an expression which describes what segment cost function
should be constructed for use with PELT. For parametric segment cost functions,
this is represented by the name of a distribution (as in the Distributions package)
with some parameters replaced by '?' to indicate that a parameters is changing.
The full list of available cost functions based on parametric distributions is
as follows:

* `Normal(?, σ): Normal model with changing mean and fixed standard deviation σ (see also NormalMeanSegment)
* `Normal(μ, ?)`: Normal model with fixed mean μ and changing standard deviation (see also NormalVarSegment)
* `Normal(?, ?)`: Normal model with changing mean and standard deviation (see also NormalMeanVarSegment)
* `Exponential(?)`: Exponential model with changing mean (see also ExponentialSegment)
* `Poisson(?)`: Poisson model with changing mean (see also PoissonSegment)
* `Gamma(?, β)`: Gamma model with fixed rate parameter β and changing shape parameter (see also GammaShapeSegment)
* `Gamma(α, ?)`: Gamma model with fixed shape parameter α and changing rate parameter (see also GammaRateSegment)

A nonparametric cost function is also provided and the model is represented with
the following expression:

* `Nonparametric(k)`: Nonparametric cost function with parameter k (see also NonparametricSegment)

# Example
```julia
n = 1000
λ = 100
μ, σ = Normal(0.0, 10.0), 1.0
# Samples changepoints from Normal distribution with changing mean
sample, cps = @changepoint_sampler n λ Normal(μ, σ)
# Create cost function
seg_cost = @segment_cost sample Normal(?, σ)
# Calculate changepoints using PELT and BS
pelt_cps, cost = PELT(seg_cost, n)
bs_cps = BS(seg_cost, n)
```
"""
macro segment_cost(data, dist)
    esc(cost_function(data, dist))
end

"""
Runs the PELT algorithm using a specified cost function and penalty value to find the position and number of changepoints

# Usage

1. `@PELT data changepoint_model`: Run PELT with default penalty value

2. `@PELT data changepoint_model β`: Run PELT at penalty value β

3. `@PELT data changepoint_model β₁ β₂`: Run CROPS algorithm for penalties between β₁ and β₂

# Example
```
n = 1000
λ = 100
μ, σ = Normal(0.0, 10.0), 1.0
# Samples changepoints from Normal distribution with changing mean
sample, cps = @changepoint_sampler n λ Normal(μ, σ)
# Run PELT on sample
pelt_cps, pelt_cost = @PELT sample Normal(?, σ)
```

# See also
PELT, @segment_cost
"""
macro PELT(data, dist, args...)
    cost_func = cost_function(data, dist)
    if length(args) == 0
        return esc(:(PELT($(cost_func), length($data))))
    elseif length(args) == 1
        return esc(:(PELT($(cost_func), length($data), pen=$(args[1]))))
    else
        return esc(:(CROPS($(cost_func), length($data), ($(args[1]), $(args[2])))))
    end
end

"""
# Description
Runs the Binary Segmentation algorithm using a specified cost function for a given penalty

# Usage
1. `@BS data changepoint_model`: Run binary segmentation with default penalty value

2. `@BS data changepoint_model β`: Run binary segmentation with penalty value β

# Example
```
n = 1000   # Length of time series
λ = 100    # Frequency of changepoints
α, β = Uniform(0.0, 10.0), 1.0
# Samples changepoints from Gamma distribution with changing shape
sample, cps = @changepoint_sampler n λ Gamma(α, β)
# Run binary segmentation on sample
bs_cps = @BS sample Gamma(?, β)
```

# See also
BS, @segment_cost
"""
macro BS(data, dist, args...)
    cost_func = cost_function(data, dist)
    if length(args) == 0
        return esc(:(BS($(cost_func), length($data))))
    else
        return esc(:(BS($(cost_func), length($data), pen=$(args[1]))))
    end
end


text:
module Clp

if haskey(ENV,"JULIA_CLP_LIBRARY_PATH") || VERSION < v"1.3"
    deps_file = joinpath(dirname(@__DIR__), "deps", "deps.jl")
    if isfile(deps_file)
        include(deps_file)
    else
        error("Clp not properly installed. Please run import `Pkg; Pkg.build(\"Clp\")`.")
    end
else
    import Clp_jll: libClp
end

using CEnum

include("gen/ctypes.jl")
include("gen/libclp_common.jl")
include("gen/libclp_api.jl")

const _CLP_VERSION = VersionNumber(
    "$(Clp_VersionMajor()).$(Clp_VersionMinor()).$(Clp_VersionRelease())"
)

if !(v"1.17.2" <= _CLP_VERSION <= v"1.17.6")
    error(
        "You have installed version $_CLP_VERSION of Clp, which is not " *
        "supported by Clp.jl. If the version change was breaking, changes " *
        "will need to be made to the Julia code. Please open an issue at " *
        "https://github.com/jump-dev/Clp.jl."
    )
end

include("MOI_wrapper/MOI_wrapper.jl")

# TODO(odow): remove at Clp.jl v1.0.0.
function ClpSolver(args...; kwargs...)
    error(
        "`ClpSolver` is no longer supported. If you are using JuMP, upgrade " *
        "to the latest version and use `Clp.Optimizer` instead. If you are " *
        "using MathProgBase (e.g., via `lingprog`), you will need to upgrade " *
        "to MathOptInterface (https://github.com/jump-dev/MathOptInterface.jl)."
    )
end
export ClpSolver

end


text:
using Oceananigans.Grids: total_extent

#####
##### Grid utilities and such
#####

function test_xnode_ynode_znode_are_correct(FT, N=3)

    grid = RegularCartesianGrid(FT, size=(N, N, N), x=(0, π), y=(0, π), z=(0, π),
                                topology=(Periodic, Periodic, Bounded))

    @test xnode(Cell, 2, grid) ≈ FT(π/2)
    @test ynode(Cell, 2, grid) ≈ FT(π/2)
    @test znode(Cell, 2, grid) ≈ FT(π/2)

    @test xnode(Face, 2, grid) ≈ FT(π/3)
    @test ynode(Face, 2, grid) ≈ FT(π/3)
    @test znode(Face, 2, grid) ≈ FT(π/3)

    @test xC(2, grid) == xnode(Cell, 2, grid)
    @test yC(2, grid) == ynode(Cell, 2, grid)
    @test zC(2, grid) == znode(Cell, 2, grid)

    @test xF(2, grid) == xnode(Face, 2, grid)
    @test yF(2, grid) == ynode(Face, 2, grid)
    @test zF(2, grid) == znode(Face, 2, grid)

    return nothing
end

#####
##### Regular Cartesian grids
#####

function regular_cartesian_correct_size(FT)
    grid = RegularCartesianGrid(FT, size=(4, 6, 8), extent=(2π, 4π, 9π))

    # Checking ≈ as the grid could be storing Float32 values.
    return (grid.Nx ≈ 4  && grid.Ny ≈ 6  && grid.Nz ≈ 8 &&
            grid.Lx ≈ 2π && grid.Ly ≈ 4π && grid.Lz ≈ 9π)
end

function regular_cartesian_correct_extent(FT)
    grid = RegularCartesianGrid(FT, size=(4, 6, 8), x=(1, 2), y=(π, 3π), z=(0, 4))
    return (grid.Lx ≈ 1 && grid.Ly ≈ 2π  && grid.Lz ≈ 4)
end

function regular_cartesian_correct_coordinate_lengths(FT)
    grid = RegularCartesianGrid(FT, size=(2, 3, 4), extent=(1, 1, 1), halo=(1, 1, 1),
                                topology=(Periodic, Bounded, Flat))

    return (
            length(grid.xC) == 4 &&
            length(grid.yC) == 5 &&
            length(grid.zC) == 6 &&
            length(grid.xF) == 4 &&
            length(grid.yF) == 6 &&
            length(grid.zF) == 6
           )
end

function regular_cartesian_correct_halo_size(FT)
    grid = RegularCartesianGrid(FT, size=(4, 6, 8), extent=(2π, 4π, 9π), halo=(1, 2, 3))
    return (grid.Hx == 1  && grid.Hy == 2  && grid.Hz == 3)
end

function regular_cartesian_correct_halo_faces(FT)
    N, H, L = 4, 1, 2.0
    Δ = L / N
    grid = RegularCartesianGrid(FT, size=(N, N, N), x=(0, L), y=(0, L), z=(0, L), halo=(H, H, H))
    return grid.xF[0] == - H * Δ && grid.yF[0] == - H * Δ && grid.zF[0] == - H * Δ
end

function regular_cartesian_correct_first_cells(FT)
    N, H, L = 4, 1, 4.0
    Δ = L / N
    grid = RegularCartesianGrid(FT, size=(N, N, N), x=(0, L), y=(0, L), z=(0, L), halo=(H, H, H))
    return grid.xC[1] == Δ/2 && grid.yC[1] == Δ/2 && grid.zC[1] == Δ/2
end

function regular_cartesian_correct_end_faces(FT)
    N, L = 4, 2.0
    Δ = L / N
    grid = RegularCartesianGrid(FT, size=(N, N, N), x=(0, L), y=(0, L), z=(0, L), halo=(1, 1, 1),
                                topology=(Periodic, Bounded, Flat))
    return grid.xF[end] == L && grid.yF[end] == L + Δ && grid.zF[end] == L
end

function regular_cartesian_ranges_have_correct_length(FT)
    Nx, Ny, Nz = 8, 9, 10
    Hx, Hy, Hz = 1, 2, 1

    grid = RegularCartesianGrid(FT, size=(Nx, Ny, Nz), extent=(1, 1, 1), halo=(Hx, Hy, Hz),
                                topology=(Bounded, Bounded, Bounded))

    return (length(grid.xC) == Nx + 2Hx && length(grid.xF) == Nx + 1 + 2Hx &&
            length(grid.yC) == Ny + 2Hy && length(grid.yF) == Ny + 1 + 2Hy &&
            length(grid.zC) == Nz + 2Hz && length(grid.zF) == Nz + 1 + 2Hz)
end

# See: https://github.com/climate-machine/Oceananigans.jl/issues/480
function regular_cartesian_no_roundoff_error_in_ranges(FT)
    Nx, Ny, Nz, Hz = 1, 1, 64, 1
    grid = RegularCartesianGrid(FT, size=(Nx, Ny, Nz), extent=(1, 1, π/2), halo=(1, 1, Hz))
    return length(grid.zC) == Nz + 2Hz && length(grid.zF) == Nz + 1 + 2Hz
end

function regular_cartesian_grid_properties_are_same_type(FT)
    grid = RegularCartesianGrid(FT, size=(10, 10, 10), extent=(1, 1//7, 2π))
    return all(isa.([grid.Lx, grid.Ly, grid.Lz, grid.Δx, grid.Δy, grid.Δz], FT)) &&
           all(eltype.([grid.xF, grid.yF, grid.zF, grid.xC, grid.yC, grid.zC]) .== FT)
end

#####
##### Vertically stretched grids
#####

function correct_constant_grid_spacings(FT)
    grid = VerticallyStretchedCartesianGrid(FT, size=(16, 16, 16), x=(0, 1), y=(0, 1), zF=collect(0:16))
    return all(grid.ΔzF .== 1) && all(grid.ΔzC .== 1)
end

function correct_quadratic_grid_spacings(FT)
    Nx = Ny = Nz = 16
    grid = VerticallyStretchedCartesianGrid(FT, size=(Nx, Ny, Nz),
                                            x=(0, 1), y=(0, 1), zF=collect(0:Nz).^2)

     zF(k) = (k-1)^2
     zC(k) = (k^2 + (k-1)^2) / 2
    ΔzF(k) = k^2 - (k-1)^2
    ΔzC(k) = 2k - 2

     zF_is_correct = all(isapprox.(  grid.zF[1:Nz+1],  zF.(1:Nz+1) ))
     zC_is_correct = all(isapprox.(  grid.zC[1:Nz],    zC.(1:Nz)   ))
    ΔzF_is_correct = all(isapprox.( grid.ΔzF[1:Nz],   ΔzF.(1:Nz)   ))

    # Note that ΔzC[1, 1, 1] involves a halo point, which is not directly determined by
    # the user-supplied zF
    ΔzC_is_correct = all(isapprox.( grid.ΔzC[2:Nz-1], ΔzC.(2:Nz-1) ))

    return zF_is_correct && zC_is_correct && ΔzF_is_correct && ΔzC_is_correct
end

function correct_tanh_grid_spacings(FT)
    Nx = Ny = Nz = 16

    S = 3  # Stretching factor

    zF(k) = tanh(S * (2 * (k - 1) / Nz - 1)) / tanh(S)

    grid = VerticallyStretchedCartesianGrid(FT, size=(Nx, Ny, Nz), x=(0, 1), y=(0, 1), zF=zF)

     zC(k) = (zF(k) + zF(k+1)) / 2
    ΔzF(k) = zF(k+1) - zF(k)
    ΔzC(k) = zC(k) - zC(k-1)

     zF_is_correct = all(isapprox.(  grid.zF[1:Nz+1],  zF.(1:Nz+1) ))
     zC_is_correct = all(isapprox.(  grid.zC[1:Nz],    zC.(1:Nz)   ))
    ΔzF_is_correct = all(isapprox.( grid.ΔzF[1:Nz],   ΔzF.(1:Nz)   ))

    # See correct_quadratic_grid_spacings for an explanation of this test component
    ΔzC_is_correct = all(isapprox.( grid.ΔzC[2:Nz-1], ΔzC.(2:Nz-1) ))

   return zF_is_correct && zC_is_correct && ΔzF_is_correct && ΔzC_is_correct
end

function vertically_stretched_grid_properties_are_same_type(FT)
    Nx, Ny, Nz = 16, 16, 16
    grid = VerticallyStretchedCartesianGrid(FT, size=(16, 16, 16), x=(0,1), y=(0,1), zF=collect(0:16))
    return all(isa.([grid.Lx, grid.Ly, grid.Lz, grid.Δx, grid.Δy], FT)) &&
           all(eltype.([grid.ΔzF, grid.ΔzC, grid.xF, grid.yF, grid.zF, grid.xC, grid.yC, grid.zC]) .== FT)
end

#####
##### Test the tests
#####

@testset "Grids" begin
    @info "Testing grids..."

    @testset "Grid utils" begin
        @info "  Testing grid utilities..."
        @test total_extent(Periodic, 1, 0.2, 1.0) == 1.2
        @test total_extent(Bounded, 1, 0.2, 1.0) == 1.4
        for FT in float_types
            test_xnode_ynode_znode_are_correct(FT)
        end
    end

    @testset "Regular Cartesian grid" begin
        @info "  Testing regular Cartesian grid..."

        @testset "Grid initialization" begin
            @info "    Testing grid initialization..."

            for FT in float_types
                @test regular_cartesian_correct_size(FT)
                @test regular_cartesian_correct_extent(FT)
                @test regular_cartesian_correct_coordinate_lengths(FT)
                @test regular_cartesian_correct_halo_size(FT)
                @test regular_cartesian_correct_halo_faces(FT)
                @test regular_cartesian_correct_first_cells(FT)
                @test regular_cartesian_correct_end_faces(FT)
                @test regular_cartesian_ranges_have_correct_length(FT)
                @test regular_cartesian_no_roundoff_error_in_ranges(FT)
                @test regular_cartesian_grid_properties_are_same_type(FT)
            end
        end

        @testset "Grid dimensions" begin
            @info "    Testing grid constructor errors..."

            for FT in float_types
                @test isbitstype(typeof(RegularCartesianGrid(FT, size=(16, 16, 16), extent=(1, 1, 1))))

                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32,), extent=(1, 1, 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32, 64), extent=(1, 1, 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32, 32, 32, 16), extent=(1, 1, 1))

                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32, 32, 32.0), extent=(1, 1, 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(20.1, 32, 32), extent=(1, 1, 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32, nothing, 32), extent=(1, 1, 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32, "32", 32), extent=(1, 1, 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32, 32, 32), extent=(1, nothing, 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32, 32, 32), extent=(1, "1", 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32, 32, 32), extent=(1, 1, 1), halo=(1, 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(32, 32, 32), extent=(1, 1, 1), halo=(1.0, 1, 1))

                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), x=2)
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), y=[1, 2])
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), z=(-π, π))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), x=1, y=2, z=3)
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), x=(0, 1), y=(0, 2), z=4)
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), x=(-1//2, 1), y=(1//7, 5//7), z=("0", "1"))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), x=(-1//2, 1), y=(1//7, 5//7), z=(1, 2, 3))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), x=(1, 0), y=(1//7, 5//7), z=(1, 2))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), x=(0, 1), y=(1, 5), z=(π, -π))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), x=(0, 1), y=(1, 5), z=(π, -π))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), extent=(1, 2, 3), x=(0, 1))
                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), extent=(1, 2, 3), x=(0, 1), y=(1, 5), z=(-π, π))

                @test_throws ArgumentError RegularCartesianGrid(FT, size=(16, 16, 16), extent=(1, 1, 1), topology=(Periodic, Periodic, Flux))
            end
        end

        # Testing show function
        topo = (Periodic, Periodic, Periodic)
        grid = RegularCartesianGrid(topology=topo, size=(3, 7, 9), x=(0, 1), y=(-π, π), z=(0, 2π))
        show(grid); println();
        @test grid isa RegularCartesianGrid
    end

    @testset "Vertically stretched Cartesian grid" begin
        @info "  Testing vertically stretched Cartesian grid..."

        @testset "Grid initialization" begin
            @info "    Testing grid initialization..."

            for FT in float_types
                @test correct_constant_grid_spacings(FT)
                @test correct_quadratic_grid_spacings(FT)
                @test correct_tanh_grid_spacings(FT)
                @test vertically_stretched_grid_properties_are_same_type(FT)
            end
        end
    end
end


text:
module App

include("c.jl")

function julia_main()::Cint
    try
        real_main()
    catch
        Base.invokelatest(Base.display_error, Base.catch_stack())
        return 1
    end
    return 0
end

if abspath(PROGRAM_FILE) == @__FILE__
    real_main()
end

end # module


text:
using Clapeyron, Test

@testset "database_lookup" begin
    params1 = Clapeyron.getparams(["water", "methanol"], ["SAFT/PCSAFT"],return_sites=false)
    @test haskey(params1, "sigma")

    # The rest of the test will be conducted with a custom dataset in the test_csvs directory.
    testspecies = ["sp1", "sp2", "sp3", "sp4", "sp5"]

    filepath_normal = ["test_csvs/normal_pair_test.csv",
                       "test_csvs/normal_single_test.csv",
                       "test_csvs/normal_assoc_test.csv",
                       "test_csvs/normal_single2_test.csv",
                       "test_csvs/normal_assoc2_test.csv"]

    filepath_clashingheaders = ["test_csvs/headercollision_single_test",
                                "test_csvs/headercollision_assoc_test"]

    filepath_asymmetry = ["test_csvs/asymmetry_pair_test",
                          "test_csvs/asymmetry_assoc_test"]

    filepath_multiple_identifiers = ["test_csvs/multiple_identifiers_single_test.csv",
                                     "test_csvs/multiple_identifiers_pair_test.csv",
                                     "test_csvs/multiple_identifiers_assoc_test.csv"]
    filepath_gc = ["test_csvs/group_test.csv"]
    filepath_param_gc = ["test_csvs/group_param_test.csv"]
    # Check that it detects the right sites.
    @test Clapeyron.findsitesincsvs(testspecies, filepath_normal) == [[],
                                                                     [],
                                                                     ["e", "e2", "H"],
                                                                     ["e", "H"],
                                                                     ["e", "e2", "H"]]

    # Check that it throws an error if ignore_missing_singleparams is not set to true.
    @test_throws ErrorException Clapeyron.getparams(testspecies; userlocations=filepath_normal,return_sites = false)

    params = Clapeyron.getparams(testspecies; userlocations=filepath_normal, ignore_missing_singleparams=["emptyparam","missingparam"],return_sites = false)
    sites = Clapeyron.SiteParam(params["intparam"].components)
    #Printing: SingleParam
    @test repr(params["intparam"]) == "SingleParam{Int64}(\"intparam\")[\"sp1\",\"sp2\",\"sp3\",\"sp4\",\"sp5\"]"
    @test repr("text/plain",params["intparam"]) == "SingleParam{Int64}(\"intparam\") with 5 components:\n \"sp1\" => 6\n \"sp2\" => 2\n \"sp3\" => 7\n \"sp4\" => 4\n \"sp5\" => 5"
    #Printing: PairParam
    @test repr(params["overwriteparam"]) == "PairParam{Float64}(\"overwriteparam\")[5×5 Matrix{Float64}]"
    @test repr("text/plain",params["overwriteparam"]) == "PairParam{PairParam{Float64}}[\"sp1\", \"sp2\", \"sp3\", \"sp4\", \"sp5\"]) with values:\n5×5 Matrix{Float64}:\n 1.6  4.0  0.0  0.0  0.0\n 4.0  1.2  0.0  3.0  0.0\n 0.0  0.0  1.3  2.0  0.0\n 0.0  3.0  2.0  1.4  0.0\n 0.0  0.0  0.0  0.0  1.5"
    #Printing: AssocParam
    @test repr(params["overwriteassocparam"]) == "AssocParam{String}(\"overwriteassocparam\")[\"val1\", \"val8\", \"val5\", \"val4\", \"val7\", \"val6\", \"val3\", \"42\"]"
    @test repr("text/plain",params["overwriteassocparam"]) == "AssocParam{String}[\"sp1\", \"sp2\", \"sp3\", \"sp4\", \"sp5\"]) with values:\n(\"sp3\", \"e\") >=< (\"sp3\", \"H\"): val1\n(\"sp3\", \"e2\") >=< (\"sp3\", \"H\"): val8\n(\"sp3\", \"H\") >=< (\"sp4\", \"e\"): val5\n(\"sp3\", \"H\") >=< (\"sp4\", \"H\"): val4\n(\"sp4\", \"e\") >=< (\"sp5\", \"H\"): val7\n(\"sp4\", \"H\") >=< (\"sp5\", \"e2\"): val6\n(\"sp5\", \"e\") >=< (\"sp5\", \"e\"): val3\n(\"sp5\", \"e\") >=< (\"sp5\", \"H\"): 42\n"
    #Printing: SiteParam
    @test repr(sites) == "SiteParam[\"sp1\" => [], \"sp2\" => [], \"sp3\" => [], \"sp4\" => [], \"sp5\" => []]"
    @test repr("text/plain",sites) == "SiteParam with 5 components:\n \"sp1\": (no sites)\n \"sp2\": (no sites)\n \"sp3\": (no sites)\n \"sp4\": (no sites)\n \"sp5\": (no sites)"
    # Check that all the types are correct.
    @test typeof(params["intparam"]) <: Clapeyron.SingleParam{Int}
    @test typeof(params["doubleparam"]) <: Clapeyron.SingleParam{Float64}
    @test typeof(params["boolparam"]) <: Clapeyron.SingleParam{Bool}
    @test typeof(params["stringparam"]) <: Clapeyron.SingleParam{String}
    # If column has both strings and numbers, they should all be strings.
    @test typeof(params["mixedparam"]) <: Clapeyron.SingleParam{String}
    # Contains missing values
    @test typeof(params["missingparam"]) <: Clapeyron.SingleParam{Int}
    # All missing values
    #before returned Clapeyron.SingleParam{Any}
    #now returns Clapeyron.SingleParam{Float64}
    @test typeof(params["emptyparam"]) <: Clapeyron.SingleParam{Float64}
    # Overwrite Int with Float, and also an upgrade from single to pair param
    @test typeof(params["overwriteparam"]) <: Clapeyron.PairParam{Float64}
    #test for missingness in PairParam, i got bitten by this before
    @test params["overwriteparam"].ismissingvalues ==  Bool[0 0 1 1 1;
                                                            0 0 1 0 1;
                                                            1 1 0 0 1;
                                                            1 0 0 0 1;
                                                            1 1 1 1 0]
    # Overwrite String with Int
    @test typeof(params["overwritestringparam"]) <: Clapeyron.SingleParam{String}
    # Overwrite Int with String
    @test typeof(params["overwriteassocparam"]) <: Clapeyron.AssocParam{String}
    @test repr(params["overwriteassocparam"]) == "AssocParam{String}(\"overwriteassocparam\")[\"val1\", \"val8\", \"val5\", \"val4\", \"val7\", \"val6\", \"val3\", \"42\"]"
    @test repr("text/plain",params["overwriteassocparam"]) == "AssocParam{String}[\"sp1\", \"sp2\", \"sp3\", \"sp4\", \"sp5\"]) with values:\n(\"sp3\", \"e\") >=< (\"sp3\", \"H\"): val1\n(\"sp3\", \"e2\") >=< (\"sp3\", \"H\"): val8\n(\"sp3\", \"H\") >=< (\"sp4\", \"e\"): val5\n(\"sp3\", \"H\") >=< (\"sp4\", \"H\"): val4\n(\"sp4\", \"e\") >=< (\"sp5\", \"H\"): val7\n(\"sp4\", \"H\") >=< (\"sp5\", \"e2\"): val6\n(\"sp5\", \"e\") >=< (\"sp5\", \"e\"): val3\n(\"sp5\", \"e\") >=< (\"sp5\", \"H\"): 42\n"
    # Check that values of "sp1" and "sp3" has been correctly overwritten.
    # "sp1" was overwritten from the same file
    # "sp3" was overwritten from a separate file, "normal_single2_test.csv"
    @test params["intparam"].values == [6, 2, 7, 4, 5]
    # Also check that "testsource1" and "testsource3" are not present, and in their place
    # we have "testsource6" and "testsource19".
    @test "testsource1" ∉ params["intparam"].sources
    @test "testsource3" ∉ params["intparam"].sources
    @test "testsource5" ∈ params["intparam"].sources
    @test "testsource19" ∈ params["intparam"].sources

    # Check that missing values have been correctly defaulted.
    @test params["missingparam"].values == [0, 2, 3, 0, 0]
    @test params["missingparam"].ismissingvalues == Bool[1, 0, 0, 1, 1]

    #test that Passing a single param to a pair param doesnt erase the missings
    singletopairparam = PairParam(params["missingparam"],"singletopairparam")
    diagmissingvalues = view(singletopairparam.ismissingvalues,1:6:25)
    @test all(diagmissingvalues .== params["missingparam"].ismissingvalues)
    diagmissingvalues .= true
    @test all(singletopairparam.ismissingvalues)

    # Check that the promotion form 1D to 2D array is succesful, with non-diagonal values present and symmetrical.
    @test params["overwriteparam"].values == [1.6  4.0  0.0  0.0  0.0
                                              4.0  1.2  0.0  3.0  0.0
                                              0.0  0.0  1.3  2.0  0.0
                                              0.0  3.0  2.0  1.4  0.0
                                              0.0  0.0  0.0  0.0  1.5]

    @test params["overwriteparam"].diagvalues == [1.6, 1.2, 1.3, 1.4, 1.5]

    assoc_param_values =
    [[Array{Int64}(undef,0,0)]  [Array{Int64}(undef,0,0)]  [Array{Int64}(undef,0,3)          ]  [Array{Int64}(undef,0,2)]  [Array{Int64}(undef,0,3)       ]
     [Array{Int64}(undef,0,0)]  [Array{Int64}(undef,0,0)]  [Array{Int64}(undef,0,3)          ]  [Array{Int64}(undef,0,2)]  [Array{Int64}(undef,0,3)       ]
     [Array{Int64}(undef,3,0)]  [Array{Int64}(undef,3,0)]  [[0 0 2000; 0 0 2700; 2000 2700 0]]  [[0 0; 0 0; 2400 2300]  ]  [[0 0 0; 0 0 0; 0 0 0]         ]
     [Array{Int64}(undef,2,0)]  [Array{Int64}(undef,2,0)]  [[0 0 2400; 0 0 2300]             ]  [[0 0; 0 0]             ]  [[0 0 2600; 0 2500 0]          ]
     [Array{Int64}(undef,3,0)]  [Array{Int64}(undef,3,0)]  [[0 0 0; 0 0 0; 0 0 0]            ]  [[0 0; 0 2500; 2600 0]  ]  [[2200 0 2100; 0 0 0; 2100 0 0]]]

    for i ∈ 1:5
        for j ∈ 1:5
            valij = assoc_param_values[i,j]
            s1,s2 = size(valij)
            testij = params["assocparam"].values[i,j]
            if iszero(s1*s2)
                continue # a compressed assoc matrix can hold more information that is discarded normally
            else
                for a ∈ 1:s1
                    for b ∈ 1:s2
                        @test valij[a,b] == testij[a,b]
                    end
                end
            end
        end
    end
  #  @test params["assocparam"].values[i,j] .== assoc_param_values[i,j] for (i,j) in zip()

    # Clashing headers between association and non-association parameters are not allowed
    @test_throws ErrorException Clapeyron.getparams(testspecies; userlocations=filepath_clashingheaders)

    # If parameter is not tagged as asymmetrical, having non-missing values across the diagonal will throw an error
    @test_throws ErrorException Clapeyron.getparams(testspecies; userlocations=filepath_asymmetry, ignore_missing_singleparams=["asymmetricpair"],return_sites = false)

    asymmetricparams = Clapeyron.getparams(testspecies; userlocations=filepath_asymmetry, asymmetricparams=["asymmetricpair", "asymmetricassoc"], ignore_missing_singleparams=["asymmetricpair"],return_sites = false)

    asymmetricparams["asymmetricpair"].values == [0.06  0.04  0.0  0.0   0.0
                                                  0.05  0.0   0.0  0.0   0.0
                                                  0.0   0.0   0.0  0.02  0.0
                                                  0.0   0.03  0.0  0.0   0.0
                                                  0.0   0.0   0.0  0.0   0.0]
    asymmetricparams["asymmetricassoc"].values ==
    [[Array{Int64}(undef,0,0)]  [Array{Int64}(undef,0,0)]  [Array{Int64}(undef,0,3)       ]  [Array{Int64}(undef,0,2) ]  [Array{Int64}(undef,0,3)    ]
     [Array{Int64}(undef,0,0)]  [Array{Int64}(undef,0,0)]  [Array{Int64}(undef,0,3)       ]  [Array{Int64}(undef,0,2) ]  [Array{Int64}(undef,0,3)    ]
     [Array{Int64}(undef,3,0)]  [Array{Int64}(undef,3,0)]  [[0 0 2800; 0 0 0; 2000 2700 0]]  [[0 2900; 0 0; 2400 2300]]  [[0 0 0; 0 0 0; 0 0 0]      ]
     [Array{Int64}(undef,2,0)]  [Array{Int64}(undef,2,0)]  [[0 0 0; 0 0 0]                ]  [[0 0; 0 0]              ]  [[0 0 0; 0 2500 0]          ]
     [Array{Int64}(undef,3,0)]  [Array{Int64}(undef,3,0)]  [[0 0 0; 0 0 0; 0 0 0]         ]  [[0 0; 0 0; 2600 0]      ]  [[2200 0 0; 0 0 0; 2100 0 0]]]

    # Also, since a non-missing value exists on the diagonal of "asymmetricpair",
    # and the diagonal contains missing values, it should throw an error
    # when parameter is not in ignore_missing_singleparams
    @test_throws ErrorException Clapeyron.getparams(testspecies; userlocations=filepath_asymmetry, asymmetricparams=["asymmetricpair", "asymmetricassoc"])

    # Testing for multiple identifiers
    multiple_identifiers = Clapeyron.getparams(testspecies; userlocations=filepath_multiple_identifiers, return_sites=false)
    @test multiple_identifiers["param_single"].values == [100, 200, 200, 300, 300]
    @test multiple_identifiers["param_pair"].values   == [1000  4000     0     0     0
                                                         4000  2000     0  6000     0
                                                            0     0  2000  5000     0
                                                            0  6000  5000  3000     0
                                                            0     0     0     0  3000]
    @test multiple_identifiers["param_assoc"].values[3,2][1,1] == 10000 &&
            multiple_identifiers["param_assoc"].values[5,1][1,1] == 20000

    # GC test, 3 comps, 4 groups

    components_gc = GroupParam(["test1", "test2", ("test3", ["grp1" => 2, "grp2" => 2, "grp3" => 3,"grp4" => 5])]; usergrouplocations=filepath_gc)

    #Printing: GroupParam
    @test repr(components_gc) == "GroupParam[\"test1\" => [\"grp1\" => 1, \"grp2\" => 2], \"test2\" => [\"grp2\" => 1], \"test3\" => [\"grp1\" => 2, \"grp2\" => 2, \"grp3\" => 3, \"grp4\" => 5]]"
    @test repr("text/plain",components_gc) == "GroupParam with 3 components:\n \"test1\": \"grp1\" => 1, \"grp2\" => 2\n \"test2\": \"grp2\" => 1\n \"test3\": \"grp1\" => 2, \"grp2\" => 2, \"grp3\" => 3, \"grp4\" => 5"

    @test components_gc.components == ["test1", "test2", "test3"]
    @test components_gc.groups == [["grp1","grp2"],["grp2"],["grp1","grp2","grp3","grp4"]]
    @test components_gc.n_groups == [[1,2], [1], [2,2,3,5]]

    # Check that flattening of groups is correct.
    @test components_gc.flattenedgroups == ["grp1", "grp2", "grp3","grp4"]
    @test components_gc.n_flattenedgroups == [[1,2,0,0], [0,1,0,0 ], [2,2,3,5]]
    @test components_gc.i_flattenedgroups == 1:4
    # Build param struct using the gc components above

    param_gc = getparams(components_gc; userlocations=filepath_param_gc)
    @test param_gc["param1"].values == [1, 2, 3, 4]

    #ParamTable
    file = ParamTable(:single,(species = ["sp1"],userparam = [2]))
    param_user = getparams(testspecies,userlocations = [file],ignore_missing_singleparams=["userparam"])
    @test param_user["userparam"].values[1] === 2
end


text:
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
#
# Description
# ==============================================================================
#
#   Functions related to the conversion from small Euler angles to quaternion.
#
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #

export smallangle_to_quat

"""
    smallangle_to_quat(θx::Number, θy::Number, θz::Number)

Create a quaternion from three small rotations of angles `θx`, `θy`, and `θz`
[rad] about the axes X, Y, and Z, respectively.

!!! note
    The quaternion is always normalized.

# Example

```jldoctest
julia> smallangle_to_quat(+0.01, -0.01, +0.01)
Quaternion{Float64}:
  + 0.999963 + 0.00499981⋅i - 0.00499981⋅j - 0.00499981⋅k
```
"""
function smallangle_to_quat(θx::Number, θy::Number, θz::Number)
    q0     = 1
    q1     = θx / 2
    q2     = θy / 2
    q3     = θz / 2
    norm_q = √(q0^2 + q1^2 + q2^2 + q3^2)

    return Quaternion(q0 / norm_q, q1 / norm_q, q2 / norm_q, q3 / norm_q)
end


text:
struct ProjectTOMLDict{T <: AbstractDict}
    project::T
end

struct ManifestTOMLDict{T <: AbstractDict}
    manifest::T
end

@enum ManifestFormatVersion manifest_format_1_0