text:
import OpenCL.cl.CLArray
import Base.LinAlg.BLAS: axpy!, scal!, gemm!

#### common stuff

function transval(t::Char)
    return t == 'T' ? clblasTrans : clblasNoTrans
end

########################### L1 functions ###############################

## AXPY
for (func, typ) in [(:clblasSaxpy, Float32), (:clblasDaxpy, Float64),
                    (:clblasCaxpy, CL_float2), (:clblasZaxpy, CL_double2)]

    @eval function axpy!(alpha::$typ, x::CLArray{$typ,1}, y::CLArray{$typ,1};
                         queue=cl.queue(x))
        @assert length(x) == length(y) "x and y have different sizes"
        return $func(Csize_t(length(x)), alpha,
                     pointer(x), Csize_t(0), Cint(1),
                     pointer(y), Csize_t(0), Cint(1), [queue])
    end

end

## SCAL
for (func, typ) in [(:clblasSscal, Float32), (:clblasDscal, Float64),
                    (:clblasCscal, CL_float2), (:clblasZscal, CL_double2)]

    @eval function scal!(n::Integer, DA::$typ,
                         DX::CLArray{$typ}, incx::Integer;
                         queue=cl.queue(DX))
        return $func(Csize_t(n), DA,
                     pointer(DX), Csize_t(0), Cint(incx),
                     [queue])
    end

end

########################### L2 functions ###############################

########################### L3 functions ###############################

## GEMM
for (func, typ) in [(:clblasSgemm, Float32), (:clblasDgemm, Float64),
                    (:clblasCgemm, CL_float2), (:clblasZgemm, CL_double2)]

    @eval function gemm!(tA::Char, tB::Char,
                         alpha::$typ, A::CLArray{$typ,2}, B::CLArray{$typ,2},
                         beta::$typ, C::CLArray{$typ,2};
                         queue=cl.queue(A))
        M = UInt64(size(A)[1])
        N = UInt64(size(B)[2])
        K = UInt64(size(A)[2])
        return $func(clblasColumnMajor, transval(tA), transval(tB), M, N, K,
              alpha, cl.pointer(A), UInt(0), M, cl.pointer(B), UInt(0), K,
              beta, cl.pointer(C), UInt(0), M, [queue])
    end

end


text:
using QuantumControl
using IOCapture
using Test

@testset "QuantumControl.jl" begin
    captured = IOCapture.capture() do
        QuantumControl.print_versions()
    end
    println(captured.output)
    @test occursin("QuantumControlBase", captured.output)
    @test occursin("Krotov", captured.output)
end

ext: jl


text:
module LIBAPS2

export APS2,
	connect!,
	disconnect!,
	enumerate_APS2s,
	init,
	set_logging_level,
	get_firmware_version,
	get_fpga_temperature,
	get_uptime,
	run,
	stop,
	trigger,
	set_run_mode,
	set_waveform_frequency,
	get_waveform_frequency,
	set_trigger_source,
	get_trigger_source,
	set_trigger_interval,
	get_trigger_interval,
	set_channel_offset,
	get_channel_offset,
	set_channel_scale,
	get_channel_scale,
	set_channel_enabled,
	get_channel_enabled,
	set_channel_delay,
	get_channel_delay,
	set_mixer_amplitude_imbalance,
	get_mixer_amplitude_imbalance,
	set_mixer_phase_skew,
	get_mixer_phase_skew,
	load_waveform,
	load_sequence,
	read_memory,
	write_memory

using Compat
import Base.run

type APS2
	ip_addr::IPv4
end

APS2() = APS2(IPv4(0))

__init__() = push!(Libdl.DL_LOAD_PATH, joinpath(dirname(@__FILE__), "../../build/"))

typealias APS2_STATUS Cint

const APS2_OK = 0

const EXTERNAL = 0
const INTERNAL = 1
const SOFTWARE = 2

const RUN_SEQUENCE = 0
const TRIG_WAVEFORM = 1
const CW_WAVEFORM = 2

function check_status(status::APS2_STATUS)
	if status != APS2_OK
		cmsg = ccall((:get_error_msg, "libaps2"), Ptr{UInt8}, (APS2_STATUS,), status)
		msg = unsafe_string(cmsg)
		error("APS2 library call failed with error: $(msg)")
	end
end

macro aps2_getter(funcName, dataType)
	funcProto = :($funcName(aps::APS2))
	funcNameBis = string(funcName)
	funcBody = quote
				val = Array($dataType, 1)
				status = ccall(($funcNameBis, "libaps2"), APS2_STATUS, (Ptr{UInt8}, Ptr{$dataType}), string(aps.ip_addr), val)
				check_status(status)
				return val[1]
			end
	Expr(:function, esc(funcProto), esc(funcBody))
end

macro aps2_channel_getter(funcName, dataType)
	funcProto = :($funcName(aps::APS2, chan::Int))
	funcNameBis = string(funcName)
	funcBody = quote
				val = Array($dataType, 1)
				status = ccall(($funcNameBis, "libaps2"), APS2_STATUS,
					(Ptr{UInt8}, Cint, Ptr{$dataType}),
					string(aps.ip_addr), chan-1, val)
				check_status(status)
				return val[1]
			end
	Expr(:function, esc(funcProto), esc(funcBody))
end

macro aps2_setter(funcName, dataType)
	funcProto = :($funcName(aps::APS2, val))
	funcNameBis = string(funcName)
	funcBody = quote
				status = ccall(($funcNameBis, "libaps2"), APS2_STATUS, (Ptr{UInt8}, $dataType), string(aps.ip_addr), val)
				check_status(status)
			end
	Expr(:function, esc(funcProto), esc(funcBody))
end

macro aps2_channel_setter(funcName, dataType)
	funcProto = :($funcName(aps::APS2, chan::Int, val))
	funcNameBis = string(funcName)
	funcBody = quote
				status = ccall(($funcNameBis, "libaps2"), APS2_STATUS,
					(Ptr{UInt8}, Cint, $dataType),
					string(aps.ip_addr), chan-1, val)
				check_status(status)
			end
	Expr(:function, esc(funcProto), esc(funcBody))
end

macro aps2_call(funcName, args...)
	funcProto = :($funcName(aps::APS2))
	argNames = [symbol("arg$ct") for ct in 1:length(args)]
	for (ct,arg) in enumerate(argNames)
		push!(funcProto.args, :($arg::$(args[ct]) ) )
	end
	funcNameBis = string(funcName)
	funcBody = quote
				status = ccall(($funcNameBis, "libaps2"), APS2_STATUS, (Ptr{UInt8},), string(aps.ip_addr), $(argNames...))
				check_status(status)
			end
	Expr(:function, esc(funcProto), esc(funcBody))
end

function get_numDevices()
	numDevices = Array(Cuint,1)
	status = ccall((:get_numDevices, "libaps2"), APS2_STATUS, (Ptr{Cuint},), numDevices)
	return numDevices[1]
end

function enumerate_APS2s()
	numDevices = get_numDevices()
	Cip_addrs = Array(Ptr{UInt8}, numDevices)
	status = ccall((:get_device_IPs, "libaps2"), APS2_STATUS, (Ptr{Ptr{UInt8}},), Cip_addrs)
	ip_addrs = Array(IPv4, 0)
	for ct = 1:numDevices
		push!(ip_addrs, IPv4(bytestring(Cip_addrs[ct])))
	end
	return (numDevices, ip_addrs)
end

function connect!(aps::APS2, ip_addr::IPv4)
	if aps.ip_addr !== IPv4(0)
		warn("Disconnecting from $(aps.ip_addr) before connection to $(ip_addr)")
		disconnect!(aps)
	end
	status = ccall((:connect_APS, "libaps2"), APS2_STATUS, (Ptr{UInt8},), string(ip_addr))
	check_status(status)
	aps.ip_addr = ip_addr
end

function disconnect!(aps::APS2)
	if aps.ip_addr !== IPv4(0)
		status = ccall((:disconnect_APS, "libaps2"), APS2_STATUS, (Ptr{UInt8},), string(aps.ip_addr))
		check_status(status)
		aps.ip_addr = IPv4(0)
	end
end

function init(aps::APS2, force=1)
	status = ccall((:init_APS, "libaps2"), APS2_STATUS, (Ptr{UInt8}, Cint), string(aps.ip_addr), Cint(force))
	check_status(status)
end

@enum(logLEVEL,
	logERROR,
	logWARNING,
	logINFO,
	logDEBUG,
	logDEBUG1,
	logDEBUG2,
	logDEBUG3,
	logDEBUG4)

set_logging_level(level) = ccall((:set_logging_level, "libaps2"), APS2_STATUS, (Cint,), level)

function get_firmware_version(aps::APS2)
	ver = Array{UInt32}(1)
	git_sha1 = Array{UInt32}(1)
	build_timestamp = Array{UInt32}(1)
	ver_string = Array{UInt8}(64)
	status = ccall((:get_firmware_version, "libaps2"),
		APS2_STATUS,
		(Cstring, Ptr{UInt32}, Ptr{UInt32}, Ptr{UInt32}, Ptr{UInt8}),
		string(aps.ip_addr), ver, git_sha1, build_timestamp, ver_string)
	check_status(status)
	ver_string[end] = '\0'
	return (ver[1], git_sha1[1], build_timestamp[1], unsafe_string(pointer(ver_string)))
end

@aps2_getter get_uptime Float64
@aps2_getter get_fpga_temperature Float32

@aps2_call run
@aps2_call stop

@aps2_setter set_run_mode Cint
@aps2_setter set_waveform_frequency Float32
@aps2_getter get_waveform_frequency Float32

@aps2_channel_setter set_channel_offset Float32
@aps2_channel_getter get_channel_offset Float32

@aps2_channel_setter set_channel_scale Float32
@aps2_channel_getter get_channel_scale Float32

@aps2_channel_setter set_channel_enabled Cint
@aps2_channel_getter get_channel_enabled Cint

@aps2_channel_setter set_channel_delay Cuint
@aps2_channel_getter get_channel_delay Cuint

@aps2_setter set_mixer_amplitude_imbalance Float32
@aps2_getter get_mixer_amplitude_imbalance Float32
@aps2_setter set_mixer_phase_skew Float32
@aps2_getter get_mixer_phase_skew Float32

@aps2_setter set_trigger_source Cint
@aps2_getter get_trigger_source Cint
@aps2_setter set_trigger_interval Float64
@aps2_getter get_trigger_interval Float64
@aps2_call trigger

function load_waveform{T<:Integer}(aps::APS2, chan, wf::Vector{T})
	status = ccall((:set_waveform_int, "libaps2"), APS2_STATUS, (Ptr{UInt8}, Cint, Ptr{Int16}, UInt32), string(aps.ip_addr), Cint(chan-1), map(Int16, wf), length(wf))
	check_status(status)
end

function load_waveform{T<:Real}(aps::APS2, chan, wf::Vector{T})
	status = ccall((:set_waveform_float, "libaps2"), APS2_STATUS, (Ptr{UInt8}, Cint, Ptr{Float32}, UInt32), string(aps.ip_addr), Cint(chan-1), map(Float32, wf), length(wf))
	check_status(status)
end

function load_sequence(aps::APS2, seqfile)
	status = ccall((:load_sequence_file, "libaps2"), APS2_STATUS, (Ptr{UInt8}, Ptr{UInt8}), string(aps.ip_addr), seqfile)
	check_status(status)
end

function read_memory(aps::APS2, addr, numWords)
	buf = Array(UInt32, numWords)
	ccall((:read_memory, "libaps2"), Cint, (Ptr{Cchar}, UInt32, Ptr{UInt32}, UInt32), string(aps.ip_addr), addr, buf, numWords)
	buf
end

write_memory(aps::APS2, addr, data::Vector{UInt32}) =
	ccall((:write_memory, "libaps2"), Cint, (Ptr{Cchar}, UInt32, Ptr{UInt32}, UInt32), string(aps.ip_addr), addr, data, length(data))

function run_DAC_BIST(aps::APS2, chan, test_vec)
	results = Array{UInt32}(8)
	passed = ccall((:run_DAC_BIST, "libaps2"),
									Cint,
									(Ptr{Cchar}, Cint, Ptr{Cshort}, Cuint, Ptr{Cuint}),
									string(aps.ip_addr), chan-1, test_vec, length(test_vec), results)
	passed != 0, results
end

function toggle_DAC_clock(aps::APS2, dac)
	status = ccall((:toggle_DAC_clock, "libaps2"), APS2_STATUS, (Ptr{UInt8}, Cint), string(aps.ip_addr), dac-1)
	check_status(status)
end

function create_test_waveform()
	#Create a test pattern with the follow pattens separated by 10ns of zero:
	# * a full scale ramp
	# * gaussians with 256/128/64/32/16/8 points
	# * +/- square pulses with 256/128/64/32/16/8 points

	wf = int16([-8192:8191])
	for ex = 8:-1:3
		xpts = linspace(-2,2, 2^ex)
		vertShift = exp(-(2 + xpts[2] - xpts[1])^2/2)
		gaussWF = int16((8191/(1-vertShift)) * (exp(-xpts.^2/2) - vertShift))
		wf = cat(1, wf, zeros(Int16, 12), gaussWF)
	end

	for ex = 8:-1:3
		wf = cat(1, wf, zeros(Int16, 12), fill(int16(8191), int(2^ex/2)), fill(int16(-8192), int(2^ex/2)))
	end
	return wf
end

end


text:
## Plots the CIR process examples.
#
# All data is loaded from cir_example.jld

using Winston, Color

include("utils.jl")

# settings
datafile = "cir_examples.jld"
figpath = "figs"
# n colors as columns of returned matrix
pcols(n) = [RGB(xi, xi, xi) for xi in linspace(0.0, 0.7, n)]
pwidth = 2
yrange = [0, 1.6]

## Simulations for different mu_tau
@readdata(datafile, "var_mu_tau", mu_taus, sig_tau, theta_tau, tau_traj, ts)
p = FramedPlot(title = "CIR trajectories for different \\tau  means",
    xlabel = "t", ylabel = "\\tau",
    xrange = [ts[1], ts[end]], yrange = yrange,
    aspect_ratio = 3/4)
nmu_tau = length(mu_taus)
ccols = pcols(nmu_tau)
pcurve = cell(nmu_tau);
for i = 1:nmu_tau
    pcurve[i] = Curve(ts, tau_traj[i,:], color=ccols[i], width=pwidth)
    setattr(pcurve[i], label=@sprintf("%4.2f", mu_taus[i]))
    add(p, pcurve[i])
end
add(p, Legend(0.1, 0.9, pcurve))
display(p)
println("-- Figure saved as $figpath/cir_examples_mu_tau.eps")
savefig(p, joinpath(figpath, "cir_examples_mu_tau.eps"))

## Simulations for different sig_tau
@readdata(datafile, "var_sig_tau", mu_tau, sig_taus, theta_tau, tau_traj, ts)
p = FramedPlot(title = "CIR trajectories for different \\tau  SDs",
    xlabel = "t", ylabel = "\\tau",
    xrange = [ts[1], ts[end]], yrange = yrange,
    aspect_ratio = 3/4)
nsig_tau = length(sig_taus)
ccols = pcols(nsig_tau)
pcurve = cell(nsig_tau);
for i = 1:nsig_tau
    pcurve[i] = Curve(ts, tau_traj[i,:], color=ccols[i], width=pwidth)
    setattr(pcurve[i], label=@sprintf("%4.2f", sig_taus[i]))
    add(p, pcurve[i])
end
add(p, Legend(0.1, 0.9, pcurve))
display(p)
println("-- Figure saved as $figpath/cir_examples_sig_tau.eps")
savefig(p, joinpath(figpath, "cir_examples_sig_tau.eps"))

## Simulations for different theta_tau
@readdata(datafile, "var_theta_tau", mu_tau, sig_tau, theta_taus, tau_traj, ts)
p = FramedPlot(title = "CIR trajectories for different \\tau  speeds",
    xlabel = "t", ylabel = "\\tau",
    xrange = [ts[1], ts[end]], yrange = yrange,
    aspect_ratio = 3/4)
ntheta_tau = length(theta_taus)
ccols = pcols(ntheta_tau)
pcurve = cell(ntheta_tau);
for i = 1:ntheta_tau
    pcurve[i] = Curve(ts, tau_traj[i,:], color=ccols[i], width=pwidth)
    setattr(pcurve[i], label=@sprintf("%4.2f", theta_taus[i]))
    add(p, pcurve[i])
end
add(p, Legend(0.1, 0.9, pcurve))
display(p)
println("-- Figure saved as $figpath/cir_examples_theta_tau.eps")
savefig(p, joinpath(figpath, "cir_examples_theta_tau.eps"))


text:
module TestEquality
using Test
using CategoricalArrays

@testset "isequal and hash for CategoricalPool" begin
    pool1 = CategoricalPool(["a", "b", "c"])
    pool1b = CategoricalPool(["a", "b", "c"])
    pool2 = CategoricalPool(["c", "a", "b"])
    pool3 = CategoricalPool(["a", "b", "c", "d"])
    pool4 = CategoricalPool(String[])

    @test isequal(pool1, pool1)
    @test pool1 == pool1
    @test pool1.equalto == C_NULL

    @test isequal(pool1, pool1b)
    @test pool1 == pool1b
    @test pool1.equalto == pointer_from_objref(pool1b)
    @test pool1b.equalto == pointer_from_objref(pool1)

    @test !isequal(pool1, pool2)
    @test pool1 != pool2
    @test pool1.equalto == pointer_from_objref(pool1b)
    @test pool2.equalto == C_NULL

    @test !isequal(pool1, pool3)
    @test pool1 != pool3
    @test pool1.equalto == pointer_from_objref(pool1b)
    @test pool3.equalto == C_NULL

    @test !isequal(pool1, pool4)
    @test pool1 != pool4
    @test pool1.equalto == pointer_from_objref(pool1b)
    @test pool4.equalto == C_NULL

    pool1c = copy(pool1b)
    @test pool1b == pool1c # To set equalto field
    @test isequal(pool1, pool1c)
    push!(pool1c, "d")
    @test pool1b.equalto == pointer_from_objref(pool1c)
    @test pool1c.equalto == C_NULL
    @test pool1b != pool1c
    @test !isequal(pool1b, pool1c)
    @test pool1b.equalto == pointer_from_objref(pool1c)
    @test pool1c.equalto == C_NULL

    pool1c = copy(pool1b)
    @test pool1b == pool1c # To set equalto field
    @test isequal(pool1, pool1c)
    get!(pool1c, "d")
    @test pool1b.equalto == pointer_from_objref(pool1c)
    @test pool1c.equalto == C_NULL
    @test pool1b != pool1c
    @test !isequal(pool1b, pool1c)
    @test pool1b.equalto == pointer_from_objref(pool1c)
    @test pool1c.equalto == C_NULL

    pool1c = copy(pool1b)
    @test pool1b == pool1c # To set equalto field
    @test isequal(pool1, pool1c)
    levels!(pool1c, ["a", "b", "c", "d"])
    @test pool1b.equalto == pointer_from_objref(pool1c)
    @test pool1c.equalto == C_NULL
    @test pool1b != pool1c
    @test !isequal(pool1, pool1c)
    @test pool1b.equalto == pointer_from_objref(pool1c)
    @test pool1c.equalto == C_NULL

    @test hash(pool1) == hash(pool1b)
    @test hash(pool1) != hash(pool2) != hash(pool3) != hash(pool4)
    @test hash(pool1, UInt(1234)) == hash(pool1b, UInt(1234)) != hash(pool1)
    @test hash(pool1, UInt(1234)) != hash(pool2, UInt(1234))
    push!(pool1, "d")
    @test hash(pool1) == hash(CategoricalArrays.hashlevels(levels(pool1)))
    get!(pool1, "e")
    @test hash(pool1) == hash(CategoricalArrays.hashlevels(levels(pool1)))

    poolnan1 = CategoricalPool([1.0, NaN])
    poolnan2 = CategoricalPool([1.0, NaN])
    @test poolnan1 == poolnan2
    @test isequal(poolnan1, poolnan2)

    poolzero1 = CategoricalPool([1.0, -0.0])
    poolzero2 = CategoricalPool([1.0, 0.0])
    @test poolzero1 != poolzero2
    @test !isequal(poolzero1, poolzero2)
end

@testset "== and isequal() for CategoricalPool{Int} and CategoricalPool{Float64}" begin
    pool1 = CategoricalPool([1, 2, 3])
    pool2 = CategoricalPool([2.0, 1.0, 3.0])

    opool1 = CategoricalPool([1, 2, 3], true)
    opool2 = CategoricalPool([2.0, 1.0, 3.0], true)

    @test isequal(pool1, pool1) === true
    @test isequal(pool1, pool2) === false
    @test isequal(pool2, pool2) === true

    @test (pool1 == pool1) === true
    @test (pool1 == pool2) === false
    @test (pool2 == pool2) === true

    @test isequal(opool1, opool1) === true
    @test isequal(opool1, opool2) === false
    @test isequal(opool2, opool2) === true

    @test (opool1 == opool1) === true
    @test (opool1 == opool2) === false
    @test (opool2 == opool2) === true

    @test (pool1 == opool1) === true
    @test (pool2 == opool2) === true
    @test (pool1 == opool2) === false

    nv1a = CategoricalValue(pool1, 1)
    nv2a = CategoricalValue(pool2, 1)
    nv1b = CategoricalValue(pool1, 2)
    nv2b = CategoricalValue(pool2, 2)

    @test isequal(nv1a, nv1a) == true
    @test isequal(nv1a, nv2a) == false
    @test isequal(nv1a, nv1b) == false
    @test isequal(nv1a, nv2b) == true

    @test isequal(nv1b, nv1a) == false
    @test isequal(nv1b, nv2a) == true
    @test isequal(nv1b, nv1b) == true
    @test isequal(nv1b, nv2b) == false

    @test isequal(nv2a, nv1a) == false
    @test isequal(nv2a, nv2a) == true
    @test isequal(nv2a, nv1b) == true
    @test isequal(nv2a, nv2b) == false

    @test isequal(nv2b, nv1a) == true
    @test isequal(nv2b, nv2a) == false
    @test isequal(nv2b, nv1b) == false
    @test isequal(nv2b, nv2b) == true

    @test isequal(1, nv1a) == true
    @test isequal(1, nv2a) == false
    @test isequal(1, nv1b) == false
    @test isequal(1, nv2b) == true

    @test isequal(nv1a, 2) == false
    @test isequal(nv2a, 2) == true
    @test isequal(nv1b, 2) == true
    @test isequal(nv2b, 2) == false

    ov1a = CategoricalValue(opool1, 1)
    ov2a = CategoricalValue(opool2, 1)
    ov1b = CategoricalValue(opool1, 2)
    ov2b = CategoricalValue(opool2, 2)

    @test isequal(ov1a, ov1a) == true
    @test isequal(ov1a, ov2a) == false
    @test isequal(ov1a, ov1b) == false
    @test isequal(ov1a, ov2b) == true

    @test isequal(ov1b, ov1a) == false
    @test isequal(ov1b, ov2a) == true
    @test isequal(ov1b, ov1b) == true
    @test isequal(ov1b, ov2b) == false

    @test isequal(ov2a, ov1a) == false
    @test isequal(ov2a, ov2a) == true
    @test isequal(ov2a, ov1b) == true
    @test isequal(ov2a, ov2b) == false

    @test isequal(ov2b, ov1a) == true
    @test isequal(ov2b, ov2a) == false
    @test isequal(ov2b, ov1b) == false
    @test isequal(ov2b, ov2b) == true

    @test isequal(1, ov1a) == true
    @test isequal(1, ov2a) == false
    @test isequal(1, ov1b) == false
    @test isequal(1, ov2b) == true

    @test isequal(ov1a, 2) == false
    @test isequal(ov2a, 2) == true
    @test isequal(ov1b, 2) == true
    @test isequal(ov2b, 2) == false

    @test (ov1a == ov1a) == true
    @test (ov1a == ov2a) == false
    @test (ov1a == ov1b) == false
    @test (ov1a == ov2b) == true

    @test (ov1b == ov1a) == false
    @test (ov1b == ov2a) == true
    @test (ov1b == ov1b) == true
    @test (ov1b == ov2b) == false

    @test (ov2a == ov1a) == false
    @test (ov2a == ov2a) == true
    @test (ov2a == ov1b) == true
    @test (ov2a == ov2b) == false

    @test (ov2b == ov1a) == true
    @test (ov2b == ov2a) == false
    @test (ov2b == ov1b) == false
    @test (ov2b == ov2b) == true

    @testset "ordered and non-ordered values are equal" begin
    @test (ov1a == nv1a) === true
    @test (ov1a == nv2a) === false
    @test (ov1a == nv1b) === false
    @test (ov1a == nv2b) === true

    @test (ov1b == nv1a) === false
    @test (ov1b == nv2a) === true
    @test (ov1b == nv1b) === true
    @test (ov1b == nv2b) === false

    @test (ov2a == nv1a) === false
    @test (ov2a == nv2a) === true
    @test (ov2a == nv1b) === true
    @test (ov2a == nv2b) === false

    @test (ov2b == nv1a) === true
    @test (ov2b == nv2a) === false
    @test (ov2b == nv1b) === false
    @test (ov2b == nv2b) === true
    end

    @testset "non-equality with missing" begin
    @test ismissing(nv1a == missing)
    @test ismissing(ov1a == missing)
    @test ismissing(missing == nv1a)
    @test ismissing(missing == ov1a)

    @test isequal(nv1a, missing) == false
    @test isequal(ov1a, missing) == false
    @test isequal(missing, nv1a) == false
    @test isequal(missing, ov1a) == false
    end
end

@testset "in()" begin
    pool = CategoricalPool([5, 1, 3])
    nv = CategoricalValue(pool, 2)

    @test (nv in 1:3) === true
    @test (nv in [1, 2, 3]) === true
    @test (nv in 2:3) === false
    @test (nv in [2, 3]) === false
end

end


text:
using KNITRO, MathProgBase
using Compat.Test

mutable struct Myquad <: MathProgBase.AbstractNLPEvaluator
end

MathProgBase.features_available(d::Myquad) = [:Grad,:HessVec]
function MathProgBase.initialize(d::Myquad,requested_features::Vector{Symbol})
    for feat in requested_features
        if !(feat in [:Grad,:Jac,:HessVec])
            error("Unsupported feature $feat")
        end
    end
end

MathProgBase.eval_g(d::Myquad,g,x) = nothing
MathProgBase.eval_jac_g(d::Myquad,g1,g2) = nothing
MathProgBase.jac_structure(d::Myquad) = Int64[], Int64[]

function MathProgBase.eval_f(d::Myquad, x)
    return (x[1]-1.)^2+(x[2]-2.)^2
end

function MathProgBase.eval_grad_f(d::Myquad, grad_f, x)
    grad_f[1] = 2*(x[1]-1.)
    grad_f[2] = 2*(x[2]-2.)
end

function MathProgBase.eval_hesslag(d::Myquad,h,x,σ,μ)
    h[1] = 2*σ
    h[2] = 2*σ
end

m = MathProgBase.NonlinearModel(KnitroSolver(KTR_PARAM_OUTLEV=2))
MathProgBase.loadproblem!(m, 2, 0, -3. * ones(2), 3. * ones(2), Float64[], Float64[], :Min, Myquad())
MathProgBase.setwarmstart!(m,rand(2))
MathProgBase.optimize!(m)
x = MathProgBase.getsolution(m)

# --- test optimal solutions ---
@testset "Test optimal solutions" begin
    @test isapprox(x[1],1.,atol=1.0e-5)
    @test isapprox(x[2],2.,atol=1.0e-5)
end


text:
# https://docs.aws.amazon.com/batch/latest/userguide/job_states.html
@doc """
    JobState

An enum for representing different possible AWS Batch job states.

See [docs](http://docs.aws.amazon.com/batch/latest/userguide/job_states.html) for details.
""" JobState

@enum JobState SUBMITTED PENDING RUNNABLE STARTING RUNNING SUCCEEDED FAILED

const STATE_MAP = Dict(string(s) => s for s in instances(JobState))
Base.parse(::Type{JobState}, str::AbstractString) = STATE_MAP[str]


text:
function jsexpr_joined(io::IO, xs, delim=", ")
  isempty(xs) && return
  for i = 1:length(xs)-1
    jsexpr(io, xs[i])
    print(io, delim)
  end
  jsexpr(io, xs[end])
end

jsexpr_joined(xs, delim=", ") = sprint(jsexpr_joined, xs, delim)

function call_expr(io, f, args...)
  f in [:(==), :+, :-, :*, :/, :%] &&
    return print(io, "(", jsexpr_joined(args, string(f)), ")")
  jsexpr(io, f)
  print(io, "(")
  jsexpr_joined(io, args)
  print(io, ")")
end

function func_expr(io, args, body; level = 0)
  named = isexpr(args, :call)
  named || print(io, "(")
  print(io, "function ")
  if named
    print(io, args.args[1])
    args = Expr(:tuple, args.args[2:end]...)
  end
  print(io, "(")
  isexpr(args, Symbol) ? print(io, args) : join(io, args.args, ", ")
  println(io, ") {")
  jsexpr(io, block(body), level = level+1)
  print(io, "  "^level, "}")
  named || print(io, ")")
end

function if_expr(io, xs; level = 0)
    if length(xs) >= 2  # we have an if
        print(io, "if (")
        jsexpr(io, xs[1])
        println(io, ") {")
        jsexpr(io, block(xs[2]), level = level+1)
        print(io, "  "^level, "}")
    end

    if length(xs) == 3  # Also have an else
        println(io, " else {")
        jsexpr(io, block(xs[3]), level = level+1)
        print(io, "  "^level, "}")
    end
end

function jsexpr(io::IO, x; level = 0)
  if isexpr(x, :block)
    for xe in x.args
      print(io, "  "^level)
      jsexpr(io, xe, level = level)
      println(io, ";")
    end
  elseif isexpr(x, :(=))
    isexpr(x.args[1], Symbol) && print(io, "let ")
    jsexpr(io, x.args[1])
    print(io, " = ")
    jsexpr(io, x.args[2], level = level)
  elseif isexpr(x, :.) && x.args[2] isa QuoteNode
    print(io, x.args[1], ".", x.args[2].value)
  elseif isexpr(x, :call)
    call_expr(io, x.args...)
  elseif isexpr(x, :tuple)
    print(io, "[")
    jsexpr_joined(io, x.args, ", ")
    print(io, "]")
  elseif @capture(x, global foo_ = bar_)
    jsexpr(io, foo)
    print(io, " = ")
    jsexpr(io, bar)
  elseif isexpr(x, :->) || isexpr(x, :function)
    func_expr(io, x.args..., level = level)
  elseif isexpr(x, :if)
    if_expr(io, x.args, level = level)
  elseif isexpr(x, :return)
    x.args[1] == nothing && return print(io, "return")
    print(io, "return ")
    jsexpr(io, x.args[1], level = level)
  elseif isexpr(x, :ref)
    jsexpr(io, x.args[1])
    print(io, "[")
    jsexpr(io, x.args[2])
    print(io, "]")
  else
    print(io, x)
  end
end

jsexpr(io::IO, x::Tuple; level = 0) = jsexpr(io, [x...], level = level)
jsexpr(io::IO, x::AbstractString; level = 0) = show(io, x)

jsexpr(x) = sprint(jsexpr, x)


text:
using Test
using TestSetExtensions
using LinearAlgebra
using Qaintessent

##==----------------------------------------------------------------------------------------------------------------------

@testset ExtendedTestSet "measurement operator construction" begin
    N = 5
    meas1 = MeasurementOperator.([3 * Matrix{Float64}(I, 2^N, 2^N), 2 * Matrix{Float64}(I, 2^N, 2^N)], ((1,2,3,4,5),))
    @test Qaintessent.check_commute(meas1) == true

    meas2 = MeasurementOperator.([X, Y], [(1,),(2,)])
    @test Qaintessent.check_commute(meas2) == true

    meas3 = MeasurementOperator(X, 2)
    meas4 = MeasurementOperator(Y, 4)
    @test Qaintessent.check_commute([meas3, meas4]) == true
end

@testset ExtendedTestSet "measurement operator helper function constuction" begin
    N = 2
    meas1 = MeasurementOperator.([X, Y], [(1,),(2,)])
    meas2 = mop.([X, Y], [(1,),(2,)])
    @test all(sparse_matrix.(meas1, (2,)) .≈ sparse_matrix.(meas2, (2,)))
end

@testset ExtendedTestSet "measurement operator abstract matrix exceptions" begin
    N = 2
    small_matrix = Matrix{Float64}(I, 2^(N)-1, 2^(N)-1)
    large_matrix = Matrix{Float64}(I, 2^(N)+1, 2^(N)+1)
    uneven_matrix = Matrix{Float64}(I, 2^(N), 2^(N+1))
    iwires = (3,4)
    non_hermitian_matrix = Qaintessent.sparse_matrix(circuit_gate(1, RxGate(0.2), 2))

    @test_throws ErrorException("Measurement operator must be a $(2^N) × $(2^N) matrix.") MeasurementOperator(small_matrix, iwires)
    @test_throws ErrorException("Measurement operator must be a $(2^N) × $(2^N) matrix.") MeasurementOperator(small_matrix, iwires)
    @test_throws ErrorException("Measurement operator must be a $(2^N) × $(2^N) matrix.") MeasurementOperator(large_matrix, iwires)
    @test_throws ErrorException("Measurement operator must be a $(2^N) × $(2^N) matrix.") MeasurementOperator(large_matrix, iwires)
    @test_throws ErrorException("Measurement operator must be a $(2^N) × $(2^N) matrix.") MeasurementOperator(uneven_matrix, iwires)
    @test_throws ErrorException("Measurement operator must be a $(2^N) × $(2^N) matrix.") MeasurementOperator(uneven_matrix, iwires)
    @test_throws ErrorException("Measurement operator must be Hermitian.") MeasurementOperator(non_hermitian_matrix, iwires)
    @test_throws ErrorException("Measurement operator must be Hermitian.") MeasurementOperator(non_hermitian_matrix, iwires)

    noncomm_matrix_1 = sparse_matrix(circuit_gate(1, HadamardGate(), 2))
    noncomm_matrix_2 = sparse_matrix(circuit_gate(1, X, 2))
    m = MeasurementOperator.([noncomm_matrix_1, noncomm_matrix_2], (iwires,))
    @test Qaintessent.check_commute(m) == false
end

@testset ExtendedTestSet "measurement operator abstract matrix exceptions" begin
    iwires = (1,)

    non_hermitian_matrix = RxGate(0.2)
    @test_throws ErrorException("Measurement operator must be Hermitian.") MeasurementOperator(non_hermitian_matrix, iwires)
    iwires = (1,)
    noncomm_gate_1 = HadamardGate()
    noncomm_gate_2 = X
    m = MeasurementOperator.([noncomm_gate_1, noncomm_gate_2], (iwires,))
    @test Qaintessent.check_commute(m) == false
end


text:
"""
    function checkDeclared_(abm,s)

Checks if an symbol or array of symbols has already been declared in the agent and returns an errors if duplications are found. Returns nothing otherwise.
"""
function checkDeclared_(abm::Agent,s::Symbol)

    if string(s)[end] == "_" &&
        error("Not possible to declare symbols with ending lower brace _. These symbols are reserved for the inner work of the module.")
    end

    if s in RESERVEDSYMBOLS
        error(s, " is a reserved symbol.")
    end

    if isdefined(Base.Math,s)
        error(s, " is a symbol from module Base. Avoid defining symbols of that library.")
    end

    if isdefined(Distributions,s)
        global dist = false
        try AgentBasedModels.eval(:(dist = Distributions.$s <: Distribution))
        catch
            Nothing
        end
        if dist
            error(s, " is a symbol from module Distributions. Avoid defining symbols of that library.")
        end
    end

    for k in keys(abm.declaredSymbols)
        if s in abm.declaredSymbols[k]
            error(s, " already defined in Agent Based model in ", k,".")
        end
    end

    # for k in keys(abm.declaredUpdates)
    #     if s in abm.declaredUpdates[k]
    #         error(s, " already defined in Agent Based model in ", k,".")
    #     end
    # end

    return
end

function checkDeclared_(abm::Agent,s::Array{Symbol})

    for i in s #Check double declarations
        if length(findall(s.==i))>1
            error("Variable ", i, " declared more than once.")
        end
        #Check if already declared
        checkDeclared_(abm,i)
    end

    return
end

"""
    function checkIsDeclared_(abm,s)

Checks if an symbol or array of symbols has already been declared in the agent and returns an errors if it is not found. Returns true otherwise.
"""
function checkIsDeclared_(abm::Agent,s::Symbol)

    found =  false

    for k in keys(abm.declaredSymbols)
        if s in abm.declaredSymbols[k]
            found = true
            break
        end
    end

    if !found
        error("Parameter ", s, " not found in the agent model.")
    end

end

function checkIsDeclared_(abm::Agent,s::Array{Symbol})

    for i in s #Check double declarations
        if length(findall(s.==i))>1
            error("Variable ", i, " declared more than once.")
        end
        checkIsDeclared_(abm,i)
    end

    return Nothing
end

"""
    function whereDeclared_(abm::Agent,s::Symbol)

Finds where the symbol s is declared, returns nothing if not found.
"""
function whereDeclared_(abm::Agent,s::Symbol)

    place =  (:NotDefined,:None)

    if isdefined(Base.Math,s)
        return (:Math,:None)
    elseif s in VALIDDISTRIBUTIONSCUDA
        return (:Distributions,:CUDAcompatible)
    elseif s in VALIDDISTRIBUTIONS
        return (:Distributions,:CPUcompatible)
    elseif isdefined(Main,s)
        return (:Main,:GPU)
    else
        if occursin("_2",string(s))
            ss = Meta.parse(string(s)[1:end-2])
            for k in keys(abm.declaredSymbols)
                if ss in abm.declaredSymbols[k]
                    place = (:Model,Meta.parse(string(k,"2")))
                    break
                end
            end
        elseif occursin("_1",string(s))
            ss = Meta.parse(string(s)[1:end-2])
            for k in keys(abm.declaredSymbols)
                if ss in abm.declaredSymbols[k]
                    place = (:Model,Meta.parse(string(k,"1")))
                    break
                end
            end
        elseif occursin("_i",string(s))
            ss = Meta.parse(string(s)[1:end-2])
            for k in keys(abm.declaredSymbols)
                if ss in abm.declaredSymbols[k]
                    place = (:Model,Meta.parse(string(k,"i")))
                    break
                end
            end
        elseif occursin("_j",string(s))
            ss = Meta.parse(string(s)[1:end-2])
            for k in keys(abm.declaredSymbols)
                if ss in abm.declaredSymbols[k]
                    place = (:Model,Meta.parse(string(k,"j")))
                    break
                end
            end
        elseif occursin("d_",string(s))
            ss = Meta.parse(string(s)[3:end])
            for k in keys(abm.declaredSymbols)
                if ss in abm.declaredSymbols[k]
                    place = (:Model,Meta.parse(string(k,"d")))
                    break
                end
            end
        else
            for k in keys(abm.declaredSymbols)
                if s in abm.declaredSymbols[k]
                    place = (:Model,Meta.parse(string(k)))
                    break
                end
            end
        end
    end

    return place
end

# function checkConsistency(amb::Agent,code::Expr,subsetDeclared)

#     #Check updated in the wrong place
#     updated,assigned,ref,symbol = symbols_(code)
#     use,declaration = validSymbolsAgent_(abm,subsetDeclared)

#     for i in symbol
#         #Check if symbol exists
#         if !(isdefined(Base.Math,i) || isdefined(Main,i) || isdefined(Distributions,j) || i in allSymbols) || isdefined(AgentBasedModels,i)
#             error("Symbol ", i, " declared in \n", code, " has not been declared of the scopes.")
#         end
#         #Check if symbol

#     for i in updated
#         place = whereDeclared_(abm,i)
#         if place == "Global" || place == "GlobalArray" || place === nothing
#             nothing
#         else
#             error(i, " has been updated in Global but was declared as ", place, " in code:\n", abm.declaredUpdates["UpdateGlobal"])
#         end
#     end
# end

max_stars_repo_name: gatocor/embryogenesisJulia

max_issues_repo_name: gatocor/embryogenesisJulia

max_forks_repo_name: gatocor/embryogenesisJulia


text:
crystal = Crystal(eye(2)u"nm", species=["Al", "O", "O"],
                  position=[1 1 1; 2 3 4], label=[:+, :-, :-])

@testset "> CrystalAtom, row $i" for i in 1:size(crystal, 1)
    atom = Crystals.CrystalAtoms.CrystalAtom(crystal, i)
    @test atom[:position] == crystal[i, :position]
    @test atom[:position, 1] == crystal[i, :position, 1]
    @test atom[:species] == crystal[i, :species]
    @test atom[:label] == crystal[i, :label]
    atom[:label] = Symbol("a$i")
    @test atom[:label] == crystal[i, :label] == Symbol("a$i")
    atom[:position, 1] = 2i - 15
    @test atom[:position, 1] == crystal[i, :position, 1] == 2i - 15
    @test names(atom) == names(crystal)
end

crystal = Crystal(eye(2)u"nm", species=["Al", "O", "O"],
                  position=[1 1 1; 2 3 4], label=[:+, :-, :-])
@test eltype(eachrow(crystal)) === Crystals.CrystalAtoms.CrystalAtom{typeof(crystal)}
@testset "> Iteration, row $i" for (i, atom) in enumerate(eachrow(crystal))
    @test atom[:position] == crystal[i, :position]
    @test atom[:species] == crystal[i, :species]
    @test atom[:label] == crystal[i, :label] != Symbol("a$i")
    atom[:label] = Symbol("a$i")
    @test atom[:label] == crystal[i, :label] == Symbol("a$i")
    atom[:position, 1] = 2i - 15
    @test atom[:position, 1] == crystal[i, :position, 1] == 2i - 15
end

crystal = Crystal(eye(2)u"nm", species=["Al", "O", "O"],
                  position=[1 1 1; 2 3 4], label=[:+, :-, :-])
@testset "> map" begin
    iterated = Any[]
    map(x->push!(iterated, iterated), eachatom(crystal))
    @test length(iterated) == length(crystal)
    @test size(iterated) == (length(crystal), )
end


text:
using FactCheck, Mongo

function with_test_db(f)
    client = MongoClient()
    dbname = randstring(15)

    try
        f(client, dbname)
    finally
        run_command(client, dbname, { "dropDatabase" => 1 })
    end
end

function itr_to_arr(itr)
    arr = {}
    for el in itr push!(arr, el) end
    arr
end

facts("Mongo") do

    with_test_db() do client, dbname

        context("Can insert and count documents") do
            insert(client, "$dbname.people", {"name" => "Brian", "age" => 30})
            insert(client, "$dbname.people", {"name" => "Lizzie", "age" => 30})
            @fact count(client, "$dbname.people") => 2
        end

        context("Can find_one document, delete it, then find nothing") do
            insert(client, "$dbname.animals", {"type" => "dog", "name" => "Spot"})

            query = {"type" => "dog"}

            @fact find_one(client, "$dbname.animals", query) => (d) -> begin
                (typeof(d) <: Dict
                 && d["type"] == "dog"
                 && d["name"] == "Spot"
                 && haskey(d, "_id"))
            end

            remove(client, "$dbname.animals", query)

            @fact find_one(client, "$dbname.animals", query) => nothing
        end

        context("Can find many documents") do
            insert(client, "$dbname.foobar", {"type" => "foo"})
            insert(client, "$dbname.foobar", {"type" => "bar"})

            docs = itr_to_arr(find(client, "$dbname.foobar"))

            @fact length(docs) => 2
            @fact typeof(docs[1]) => t -> t <: Dict
        end

        context("Can update documents") do
            insert(client, "$dbname.plants", {"type" => "rose", "color" => "red"})

            query = {"type" => "rose"}

            @fact find_one(client, "$dbname.plants", query) => (d) -> begin
                d["color"] == "red"
            end

            update(client, "$dbname.plants", query, {"\$set" => {"color" => "blue"}})

            @fact find_one(client, "$dbname.plants", query) => (d) -> begin
                d["color"] == "blue"
            end
        end

    end

end


text:
import Base: put!, take!, close
import DandelionWebSockets: AbstractWSClient, ClientLogicInput, send_text, stop,
                        handle, WebSocketHandler, wsconnect

export ThrottledWSClient

typealias WSCall Tuple{Function, Vector{Any}}

# ThrottledWSClient is a throttling proxy for an AbstractWSClient. It ensures that
# messages are not sent too often.
immutable ThrottledWSClient <: AbstractWSClient
    ws::AbstractWSClient
    chan::Channel{WSCall}

    function ThrottledWSClient(ws::AbstractWSClient, interval::Float64; capacity::Int=32)
        chan = Channel{WSCall}(capacity)

        @schedule begin
            for (f, args) in chan
                f(args...)
                sleep(interval)
            end
        end

        new(ws, chan)
    end
end

# Close requests are not throttled, because this doesn't send a message over to Slack.
# Note: This might lead to enqueued messages never being sent.
function stop(ws::ThrottledWSClient)
    stop(ws.ws)
    close(ws.chan)
end

send_text(ws::ThrottledWSClient, text::UTF8String) = put!(ws.chan, (send_text, [ws.ws, text]))
wsconnect(client::ThrottledWSClient, uri::Requests.URI, handler::WebSocketHandler) =
    wsconnect(client.ws, uri, handler)


text:
# This file is auto-generated by AWSMetadata.jl
using AWS
using AWS.AWSServices: amplifyuibuilder
using AWS.Compat
using AWS.UUIDs

"""
    create_component(app_id, component_to_create, environment_name)
    create_component(app_id, component_to_create, environment_name, params::Dict{String,<:Any})

Creates a new component for an Amplify app.

# Arguments
- `app_id`: The unique ID of the Amplify app to associate with the component.
- `component_to_create`: Represents the configuration of the component to create.
- `environment_name`: The name of the backend environment that is a part of the Amplify app.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"clientToken"`: The unique client token.
"""
function create_component(
    appId,
    componentToCreate,
    environmentName;
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "POST",
        "/app/$(appId)/environment/$(environmentName)/components",
        Dict{String,Any}(
            "componentToCreate" => componentToCreate, "clientToken" => string(uuid4())
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function create_component(
    appId,
    componentToCreate,
    environmentName,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "POST",
        "/app/$(appId)/environment/$(environmentName)/components",
        Dict{String,Any}(
            mergewith(
                _merge,
                Dict{String,Any}(
                    "componentToCreate" => componentToCreate,
                    "clientToken" => string(uuid4()),
                ),
                params,
            ),
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    create_theme(app_id, environment_name, theme_to_create)
    create_theme(app_id, environment_name, theme_to_create, params::Dict{String,<:Any})

Creates a theme to apply to the components in an Amplify app.

# Arguments
- `app_id`: The unique ID of the Amplify app associated with the theme.
- `environment_name`: The name of the backend environment that is a part of the Amplify app.
- `theme_to_create`: Represents the configuration of the theme to create.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"clientToken"`: The unique client token.
"""
function create_theme(
    appId, environmentName, themeToCreate; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "POST",
        "/app/$(appId)/environment/$(environmentName)/themes",
        Dict{String,Any}(
            "themeToCreate" => themeToCreate, "clientToken" => string(uuid4())
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function create_theme(
    appId,
    environmentName,
    themeToCreate,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "POST",
        "/app/$(appId)/environment/$(environmentName)/themes",
        Dict{String,Any}(
            mergewith(
                _merge,
                Dict{String,Any}(
                    "themeToCreate" => themeToCreate, "clientToken" => string(uuid4())
                ),
                params,
            ),
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    delete_component(app_id, environment_name, id)
    delete_component(app_id, environment_name, id, params::Dict{String,<:Any})

Deletes a component from an Amplify app.

# Arguments
- `app_id`: The unique ID of the Amplify app associated with the component to delete.
- `environment_name`: The name of the backend environment that is a part of the Amplify app.
- `id`: The unique ID of the component to delete.

"""
function delete_component(
    appId, environmentName, id; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "DELETE",
        "/app/$(appId)/environment/$(environmentName)/components/$(id)";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function delete_component(
    appId,
    environmentName,
    id,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "DELETE",
        "/app/$(appId)/environment/$(environmentName)/components/$(id)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    delete_theme(app_id, environment_name, id)
    delete_theme(app_id, environment_name, id, params::Dict{String,<:Any})

Deletes a theme from an Amplify app.

# Arguments
- `app_id`: The unique ID of the Amplify app associated with the theme to delete.
- `environment_name`: The name of the backend environment that is a part of the Amplify app.
- `id`: The unique ID of the theme to delete.

"""
function delete_theme(
    appId, environmentName, id; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "DELETE",
        "/app/$(appId)/environment/$(environmentName)/themes/$(id)";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function delete_theme(
    appId,
    environmentName,
    id,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "DELETE",
        "/app/$(appId)/environment/$(environmentName)/themes/$(id)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    exchange_code_for_token(provider, request)
    exchange_code_for_token(provider, request, params::Dict{String,<:Any})

Exchanges an access code for a token.

# Arguments
- `provider`: The third-party provider for the token. The only valid value is figma.
- `request`: Describes the configuration of the request.

"""
function exchange_code_for_token(
    provider, request; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "POST",
        "/tokens/$(provider)",
        Dict{String,Any}("request" => request);
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function exchange_code_for_token(
    provider,
    request,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "POST",
        "/tokens/$(provider)",
        Dict{String,Any}(mergewith(_merge, Dict{String,Any}("request" => request), params));
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    export_components(app_id, environment_name)
    export_components(app_id, environment_name, params::Dict{String,<:Any})

Exports component configurations to code that is ready to integrate into an Amplify app.

# Arguments
- `app_id`: The unique ID of the Amplify app to export components to.
- `environment_name`: The name of the backend environment that is a part of the Amplify app.

"""
function export_components(
    appId, environmentName; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "GET",
        "/export/app/$(appId)/environment/$(environmentName)/components";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function export_components(
    appId,
    environmentName,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "GET",
        "/export/app/$(appId)/environment/$(environmentName)/components",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    export_themes(app_id, environment_name)
    export_themes(app_id, environment_name, params::Dict{String,<:Any})

Exports theme configurations to code that is ready to integrate into an Amplify app.

# Arguments
- `app_id`: The unique ID of the Amplify app to export the themes to.
- `environment_name`: The name of the backend environment that is part of the Amplify app.

"""
function export_themes(
    appId, environmentName; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "GET",
        "/export/app/$(appId)/environment/$(environmentName)/themes";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function export_themes(
    appId,
    environmentName,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "GET",
        "/export/app/$(appId)/environment/$(environmentName)/themes",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    get_component(app_id, environment_name, id)
    get_component(app_id, environment_name, id, params::Dict{String,<:Any})

Returns an existing component for an Amplify app.

# Arguments
- `app_id`: The unique ID of the Amplify app.
- `environment_name`: The name of the backend environment that is part of the Amplify app.
- `id`: The unique ID of the component.

"""
function get_component(
    appId, environmentName, id; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "GET",
        "/app/$(appId)/environment/$(environmentName)/components/$(id)";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function get_component(
    appId,
    environmentName,
    id,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "GET",
        "/app/$(appId)/environment/$(environmentName)/components/$(id)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    get_theme(app_id, environment_name, id)
    get_theme(app_id, environment_name, id, params::Dict{String,<:Any})

Returns an existing theme for an Amplify app.

# Arguments
- `app_id`: The unique ID of the Amplify app.
- `environment_name`: The name of the backend environment that is part of the Amplify app.
- `id`: The unique ID for the theme.

"""
function get_theme(
    appId, environmentName, id; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "GET",
        "/app/$(appId)/environment/$(environmentName)/themes/$(id)";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function get_theme(
    appId,
    environmentName,
    id,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "GET",
        "/app/$(appId)/environment/$(environmentName)/themes/$(id)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    list_components(app_id, environment_name)
    list_components(app_id, environment_name, params::Dict{String,<:Any})

Retrieves a list of components for a specified Amplify app and backend environment.

# Arguments
- `app_id`: The unique ID for the Amplify app.
- `environment_name`: The name of the backend environment that is a part of the Amplify app.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"maxResults"`: The maximum number of components to retrieve.
- `"nextToken"`: The token to request the next page of results.
"""
function list_components(
    appId, environmentName; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "GET",
        "/app/$(appId)/environment/$(environmentName)/components";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function list_components(
    appId,
    environmentName,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "GET",
        "/app/$(appId)/environment/$(environmentName)/components",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    list_themes(app_id, environment_name)
    list_themes(app_id, environment_name, params::Dict{String,<:Any})

Retrieves a list of themes for a specified Amplify app and backend environment.

# Arguments
- `app_id`: The unique ID for the Amplify app.
- `environment_name`: The name of the backend environment that is a part of the Amplify app.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"maxResults"`: The maximum number of theme results to return in the response.
- `"nextToken"`: The token to request the next page of results.
"""
function list_themes(
    appId, environmentName; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "GET",
        "/app/$(appId)/environment/$(environmentName)/themes";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function list_themes(
    appId,
    environmentName,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "GET",
        "/app/$(appId)/environment/$(environmentName)/themes",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    refresh_token(provider, refresh_token_body)
    refresh_token(provider, refresh_token_body, params::Dict{String,<:Any})

Refreshes a previously issued access token that might have expired.

# Arguments
- `provider`: The third-party provider for the token. The only valid value is figma.
- `refresh_token_body`: Information about the refresh token request.

"""
function refresh_token(
    provider, refreshTokenBody; aws_config::AbstractAWSConfig=global_aws_config()
)
    return amplifyuibuilder(
        "POST",
        "/tokens/$(provider)/refresh",
        Dict{String,Any}("refreshTokenBody" => refreshTokenBody);
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function refresh_token(
    provider,
    refreshTokenBody,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "POST",
        "/tokens/$(provider)/refresh",
        Dict{String,Any}(
            mergewith(
                _merge, Dict{String,Any}("refreshTokenBody" => refreshTokenBody), params
            ),
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    update_component(app_id, environment_name, id, updated_component)
    update_component(app_id, environment_name, id, updated_component, params::Dict{String,<:Any})

Updates an existing component.

# Arguments
- `app_id`: The unique ID for the Amplify app.
- `environment_name`: The name of the backend environment that is part of the Amplify app.
- `id`: The unique ID for the component.
- `updated_component`: The configuration of the updated component.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"clientToken"`: The unique client token.
"""
function update_component(
    appId,
    environmentName,
    id,
    updatedComponent;
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "PATCH",
        "/app/$(appId)/environment/$(environmentName)/components/$(id)",
        Dict{String,Any}(
            "updatedComponent" => updatedComponent, "clientToken" => string(uuid4())
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function update_component(
    appId,
    environmentName,
    id,
    updatedComponent,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "PATCH",
        "/app/$(appId)/environment/$(environmentName)/components/$(id)",
        Dict{String,Any}(
            mergewith(
                _merge,
                Dict{String,Any}(
                    "updatedComponent" => updatedComponent, "clientToken" => string(uuid4())
                ),
                params,
            ),
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    update_theme(app_id, environment_name, id, updated_theme)
    update_theme(app_id, environment_name, id, updated_theme, params::Dict{String,<:Any})

Updates an existing theme.

# Arguments
- `app_id`: The unique ID for the Amplify app.
- `environment_name`: The name of the backend environment that is part of the Amplify app.
- `id`: The unique ID for the theme.
- `updated_theme`: The configuration of the updated theme.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"clientToken"`: The unique client token.
"""
function update_theme(
    appId,
    environmentName,
    id,
    updatedTheme;
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "PATCH",
        "/app/$(appId)/environment/$(environmentName)/themes/$(id)",
        Dict{String,Any}("updatedTheme" => updatedTheme, "clientToken" => string(uuid4()));
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function update_theme(
    appId,
    environmentName,
    id,
    updatedTheme,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return amplifyuibuilder(
        "PATCH",
        "/app/$(appId)/environment/$(environmentName)/themes/$(id)",
        Dict{String,Any}(
            mergewith(
                _merge,
                Dict{String,Any}(
                    "updatedTheme" => updatedTheme, "clientToken" => string(uuid4())
                ),
                params,
            ),
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end


text:
using DataFrames, CSV, Statistics, StatsBase

df = CSV.read("./preprocessing/input/land-pressure.csv")
vars = CSV.read("./preprocessing/input/land-pressure-attributes.csv")

# remove last (empty) column
select!(df, Not(names(df)[end]))

# transform wide → long and remove rows without value
df = stack(df, 6:size(df)[2])
df.year = parse.(Int, string.(df.variable))
select!(df, Not(:variable))
rename!(df, [:model, :scenario, :region, :variable, :unit, :value, :year])
df = df[.!ismissing.(df.value),:]

# limit years between 1960 and 2019
df = df[.&(df.year .>= 1960, df.year .<= 2019), :]

df = unstack(df, :year, :value)

# replace missing values with empty strings
df[ismissing.(df.scenario), :scenario] .= ""
df[ismissing.(df.model), :model] .= ""

# calc afolu emissions
emissionsAgr = stack(df[occursin.(r"\+\|Agriculture", df.variable), :], 6:size(df)[2]; variable_name = :year)
emissionsLand = stack(df[occursin.(r"\+\|Land-use Change", df.variable), :], 6:size(df)[2]; variable_name = :year)

emissionsAgr.value += emissionsLand.value
agr = unstack(emissionsAgr, :year, :value)
df = df[.!occursin.(r"Emissions", df.variable), :]
df = vcat(agr, df)

# fix fish prices
foreach(v -> df[df.variable .== "Prices|Agriculture|Fish", Symbol(v)] .= missing, 1960:1979)
df
# sort based on var-selection
df = join(vars,df,on = :variable, kind= :outer)
df = df[.!occursin.(r"Land-use Change", df.variable), :]
CSV.write("./src/assets/data/land-pressure.csv", df)

max_stars_repo_name: SensesProject/land-use

max_issues_repo_name: SensesProject/land-use

max_forks_repo_name: SensesProject/land-use


text:
export  convert2affinePCE,
        calculateAffinePCE,
        assign2multi,
        evaluatePCE,
        sampleMeasure,
        samplePCE,
        mean,
        var,
        std

"""
    calculateAffinePCE(α::Vector{Float64})::Vector{Float64}
Computes the affine PCE coefficients ``x_0`` and ``x_1`` from recurrence coefficients ``\alpha``.
"""
function calculateAffinePCE(α::Vector{Float64})
    @assert length(α)>=1 "not enough recursion coefficients"
    return [α[1], 1.]
end

calculateAffinePCE(op::OrthoPoly) = calculateAffinePCE(op.α)
calculateAffinePCE(opq::OrthoPolyQ) = calculateAffinePCE(opq.op)

function calculateAffinePCE(i::Int64,mop::MultiOrthoPoly)
    p = length(mop.name)
    @assert i<=p "mop is $p-variate, but the PC for the $i-variate entry was requested"
    calculateAffinePCE(mop.uni[i])
end

"""
```
convert2affinePCE(a::Vector{Float64},α0::Float64)
convert2affinePCE(name::String,p1::Float64,p2::Float64,d::Dict=Dict();kind::Symbol=:lbub)
```
Computes the affine PCE coefficients ``x_0`` and ``x_1`` from
```math
X = a_1 + a_2 \\Xi = x_0 + x_1 \\phi_1(\\Xi),
```
where ``\\phi_1(t) = t-\\alpha_0`` is the first-order monic basis polynomial.

For classical polynomials the `name` can be given directly. The keyword `kind in [:lbub, :μσ]`
specifies whether `p1` and `p2` have the meaning of lower/upper bounds or mean/standard deviation.
"""
function convert2affinePCE(a1::Float64,a2::Float64,α0::Float64)
    [
        a1+α0*a2;
        a2
    ]
end
function convert2affinePCE(name::String,p1::Float64,p2::Float64,d::Dict=Dict();kind::Symbol=:lbub)
    @assert kind in [:lbub, :μσ]
    name = lowercase(name)
    a1, a2, α = 0., 0., zeros(Float64,1)
    name=="beta01" && d[:shape_a]==d[:shape_b]==1 ? name="uniform01" : ()

    if name=="gaussian"
        # p1 --> mean
        # p2 --> standard deviation
        @assert p2>=0. "standard deviation has to be non-negative"
        abs(p2)<1e-4 ? (@warn "standard deviation close to zero: $p2") : ()
        a1, a2 = p1, p2
        α,~ = rm_hermite_prob(1)
    elseif name=="uniform01"
        α,~ = rm_legendre01(1)
        if kind==:lbub
            # p1 --> lower bound
            # p2 --> upper bound
            @assert p1<p2 "inconsistent bounds"
            a1, a2 = p1, p2-p1
        elseif kind==:μσ
            # p1 --> mean
            # p2 --> standard deviation
            @assert p2>=0. "standard deviation has to be nonnegative"
            abs(p2)<1e-4 ? (@warn "standard deviation close to zero: $p2") : ()
            a1, a2 = p1-sqrt(3)*p2 , 2*sqrt(3)*p2
        else
            error("not implemented")
        end
    elseif name=="beta01"
        s_α, s_β = d[:shape_a], d[:shape_b]
        α,~ = rm_jacobi01(1,Float64(s_β)-1,Float64(s_α)-1)
        if kind==:lbub
            # p1 --> lower bound
            # p2 --> upper bound
            @assert p1<p2 "inconsistent bounds"
            a1, a2 = p1, p2-p1
        elseif kind==:μσ
            # p1 --> mean
            # p2 --> standard deviation
            @assert p2>=0. "standard deviation has to be nonnegative"
            abs(p2)<1e-4 ? (@warn "standard deviation close to zero: $p2") : ()
            a1 = p1-sqrt(s_α/s_β)*sqrt(1+s_α+s_β)*p2
            a2 = (s_α+s_β)*sqrt((s_α+s_β+1)/(s_α*s_β))*p2
        else
            error("not implemented")
        end
    elseif name=="gamma"
        error("$name is not implemented")
    elseif name=="logistic"
        @assert p2>0
        α,~ = rm_logistic(1)
        a1, a2 = p1, p2
    else
        error("$name is not implemented")
    end
    convert2affinePCE(a1,a2,α[1])
end
# convert2affinePCE(name::String,p1::Float64,p2::Float64;kind::Symbol=:lbub) = convert2affinePCE(name,p1,p2,Dict();kind=kind)

function assign2multi(x::Vector{Float64},i::Int64,ind::Matrix{Int64})
    l,p = size(ind)
    nx, deg = length(x), ind[end,end]
    @assert nx<=deg+1 "inconsistent number of coefficients ($nx vs $(deg+1))"
    @assert i<=p "basis is $p-variate, you requested $i-variate"
    myind::Vector{Int64} = findUnivariateIndices(i,ind)[1:nx]
    y = spzeros(Float64,l)
    y[myind] = x
    return y
end

function sampleMeasure_byName(n::Int64,name::String,d::Dict=Dict())::Vector{Float64}
    @assert n >= 1 "invalid number $n of samples."
    name = lowercase(name)
    name == "gaussian" && return rand(Normal(),n)
    name == "uniform01" && return rand(Uniform(),n)
    name == "beta01" && return rand(Beta(d[:shape_a],d[:shape_b]),n)
    name == "gamma" && return rand(Gamma(d[:shape],1/d[:rate]),n)
    name == "logistic" && return rand(Logistic(),n)
    throw(error("$name is not implemented"))
end

function sampleMeasure_byFun(n::Int64,w::Function,dom::Tuple{Float64,Float64};method::String="adaptiverejection")
    @assert n>=1 "invalid number $n of samples."
    method = lowercase(method)
    @assert method in [ "adaptiverejection", "inversecdf"] "method $method not implemented"
    if method == "adaptiverejection"
        # works only if w is log-concave
        sampler = RejectionSampler(w,dom)
        return run_sampler!(sampler,n)
    elseif method == "rejection"
        # all purpose method but needs a solid envelope PDF
        error("method $method not yet implemented")
    elseif method == "inversecdf"
        error("method $method not yet implemented")
        # requires CDF and reliable root-finding method
    end
end

"""
__Univariate__
```
sampleMeasure(n::Int64,name::String,w::Function,dom::Tuple{Float64,Float64},symm::Bool,d::Dict;method::String="adaptiverejection")
sampleMeasure(n::Int64,m::Measure;method::String="adaptiverejection")
sampleMeasure(n::Int64,op::OrthoPoly;method::String="adaptiverejection")
sampleMeasure(n::Int64,opq::OrthoPolyQ;method::String="adaptiverejection")
```
Draw `n` samples from the measure `m` described by its
- `name`
- weight function `w`,
- domain `dom`,
- symmetry property `symm`,
- and, if applicable, parameters stored in the dictionary `d`.
By default an adaptive rejection sampling method is used (from [AdaptiveRejectionSampling.jl](https://github.com/mauriciogtec/AdaptiveRejectionSampling.jl)),
unless it is a common random variable for which [Distributions.jl](https://github.com/JuliaStats/Distributions.jl) is used.

The function is multiply dispatched to accept `OrthoPoly` or `OrthoPolyQ`.

__Multivariate__
```
sampleMeasure(n::Int64,m::MultiMeasure;method::Vector{String}=["adaptiverejection" for i=1:length(m.name)])
sampleMeasure(n::Int64,mop::MultiOrthoPoly;method::Vector{String}=["adaptiverejection" for i=1:length(mop.meas.name)])
```
Multivariate extension which provides array of samples with `n` rows and
as many columns as the multimeasure has univariate measures.
"""
function sampleMeasure(n::Int64,name::String,w::Function,dom::Tuple{Float64,Float64},symm::Bool,d::Dict;method::String="adaptiverejection")
    # ToDo: symmetry is currently not exploited, but could be for rejection sampling
    name = lowercase(name)
    if name in ["gaussian","uniform01","beta01","gamma","logistic"]
        sampleMeasure_byName(n,name,d)
    else
        sampleMeasure_byFun(n,w,dom;method=method)
    end
end
sampleMeasure(n::Int64,m::Measure;method::String="adaptiverejection") = sampleMeasure(n,m.name,m.w,m.dom,m.symmetric,m.pars;method=method)
sampleMeasure(n::Int64,op::OrthoPoly;method::String="adaptiverejection") = sampleMeasure(n,op.meas;method=method)
sampleMeasure(n::Int64,opq::OrthoPolyQ;method::String="adaptiverejection") = sampleMeasure(n::Int64,opq.op;method=method)

function sampleMeasure(n::Int64, name::Vector{String}, w::Vector{Function},
                       dom::Vector{Tuple{Float64,Float64}}, symm::BitArray,d::Vector{Dict};
                       method::Vector{String}=["adaptiverejection" for i=1:length(name)])
    @assert length(name)==length(w)==length(dom)==length(symm)==length(d) "inconsistent number of parameters"
    m = length(name)
    ξ = zeros(n,m)
    for k=1:m
        ξ[:,k] = sampleMeasure(n,name[k],w[k],dom[k],symm[k],d[k];method=method[k])
    end
    return ξ
end
function sampleMeasure(n::Int64,m::MultiMeasure;method::Vector{String}=["adaptiverejection" for i=1:length(m.name)])
    s = falses(length(m.name))
    s[:] = m.symmetric[:]
    sampleMeasure(n,m.name,m.w_uni,m.dom,s,m.pars;method=method)
end
sampleMeasure(n::Int64,mop::MultiOrthoPoly;method::Vector{String}=["adaptiverejection" for i=1:length(mop.meas.name)]) = sampleMeasure(n,mop.meas;method=method)

"""
    evaluatePCE(x::Vector{Float64},ξ::Vector{Float64},α::Vector{Float64},β::Vector{Float64})
Evaluation of polynomial chaos expansion
```math
\\mathsf{x} = \\sum_{i=0}^{L} x_i \\phi_i{\\xi_j},
```
where `L+1 = length(x)` and ``x_j`` is the ``j``th sample where ``j=1,\\dots,m``
with `m = length(ξ)`.
"""
function evaluatePCE(x::Vector{Float64},ξ::Vector{Float64},α::Vector{Float64},β::Vector{Float64})
    @assert length(α) == length(β) "inconsistent number of recurrence coefficients"
    Nsmpl = length(ξ)
    @assert Nsmpl>=1 "inconsistent number $Nsmpl of samples"
    Nx, Nrec = length(x), length(α)
    @assert Nx<=Nrec "not enough recursion coefficients"
    if Nrec > Nx α,β = α[1:Nx], β[1:Nx] end
    ϕ = zeros(Float64,Nsmpl,Nx)
    for n in 1:Nx
        ϕ[:,n] = evaluate(n-1,ξ,α,β)
    end
    ϕ*x
end
evaluatePCE(x::Vector{Float64},ξ::Float64,α::Vector{Float64},β::Vector{Float64}) = evaluatePCE(x,[ξ],α,β)
evaluatePCE(x::Vector{Float64},ξ::Vector{Float64},op::OrthoPoly) = evaluatePCE(x,ξ,op.α,op.β)
evaluatePCE(x::Vector{Float64},ξ::Vector{Float64},opq::OrthoPolyQ) = evaluatePCE(x,ξ,opq.op)

function evaluatePCE(x::Vector{Float64},ξ::Matrix{Float64},α::Vector{Vector{Float64}},β::Vector{Vector{Float64}},ind::Matrix{Int64})
    Nsmpl = size(ξ,1)
    @assert Nsmpl>=1 "inconsistent number $Nsmpl of samples"
    @assert length(α)==length(β)==size(ξ,2)==size(ind,2) "inconsistent number of coefficients"
    Nx = length(x)
    @assert Nx <= size(ind,1) "too few pc coefficients (resp: too small basis)"
    ϕ = zeros(Float64,Nsmpl,Nx)
    for n in 1:Nx
        ϕ[:,n] = evaluate(ind[n,:],ξ,α,β)
    end
    ϕ*x
end
function evaluatePCE(x::Vector{Float64},ξ::Matrix{Float64},mOP::MultiOrthoPoly)
    a,b = coeffs(mOP)
    evaluatePCE(x,ξ,a,b,mOP.ind)
end

"""
__Univariate__
```
samplePCE(n::Int64,x::Vector{Float64},op::OrthoPoly;method::String="adaptiverejection")
samplePCE(n::Int64,x::Vector{Float64},opq::OrthoPolyQ;method::String="adaptiverejection")
```
Combines [`sampleMeasure`](@ref) and [`evaluatePCE`](@ref), i.e. it first draws `n` samples
from the measure, then evaluates the PCE for those samples.

__Multivariate__
```
samplePCE(n::Int64,x::Vector{Float64},mop::MultiOrthoPoly;method::Vector{String}=["adaptiverejection" for i=1:length(mop.meas.name)])
```
"""
function samplePCE(n::Int64,x::Vector{Float64},op::OrthoPoly;method::String="adaptiverejection")
    ξ = sampleMeasure(n,op;method=method)
    evaluatePCE(x,ξ,op)
end

samplePCE(n::Int64,x::Vector{Float64},opq::OrthoPolyQ;method::String="adaptiverejection") = samplePCE(n,x,opq.op;method=method)

function samplePCE(n::Int64,x::Vector{Float64},mop::MultiOrthoPoly;method::Vector{String}=["adaptiverejection" for i=1:length(mop.meas.name)])
    ξ = sampleMeasure(n,mop;method=method)
    evaluatePCE(x,ξ,mop)
end

"""
__Univariate__
```
mean(x::Vector{Float64},op::OrthoPoly)
mean(x::Vector{Float64},opq::OrthoPolyQ)
```
__Multivariate__
```
mean(x::Vector{Float64},mop::MultiOrthoPoly)
```
compute mean of random variable with PCE `x`
"""
function mean(x::Vector{Float64},op::OrthoPoly)::Float64
    x[1]*computeSP2(0,op.β)
end
mean(x::Vector{Float64},opq::OrthoPolyQ)::Float64 = mean(x,opq.op)
function mean(x::Vector{Float64},mop::MultiOrthoPoly)::Float64
    nunc = length(mop.uni)
    x[1]*computeSP(zeros(Int64,nunc),mop)
end

"""
__Univariate__
```
var(x::Vector{Float64},op::OrthoPoly)
var(x::Vector{Float64},opq::OrthoPolyQ)
```
__Multivariate__
```
var(x::Vector{Float64},mop::MultiOrthoPoly)
```
compute variance of random variable with PCE `x`
"""
function var(x::Vector{Float64},op::OrthoPoly)::Float64
    t2 = computeSP2(op)
    @assert length(t2)>=length(x) "cannot compute variance; too many PCE coefficients"
    sum( x[i]^2*t2[i] for i=2:length(x) )
end
var(x::Vector{Float64},opq::OrthoPolyQ)::Float64 = var(x,opq.op)

function var(x::Vector{Float64},mop::MultiOrthoPoly)::Float64
    L = size(mop.ind,1)
    @assert length(x)<=L "cannot compute variance; too many PCE coefficients"
    t2 = Tensor(2,mop)
    sum( x[i]^2*t2.get([i-1,i-1]) for i=2:length(x))
end

"""
__Univariate__
```
std(x::Vector{Float64},op::OrthoPoly)
std(x::Vector{Float64},opq::OrthoPolyQ)
```
__Multivariate__
```
std(x::Vector{Float64},mop::MultiOrthoPoly)
```
compute standard deviation of random variable with PCE `x`
"""
function std(x::Vector{Float64},op::OrthoPoly)::Float64
    sqrt(var(x,op))
end
std(x::Vector{Float64},opq::OrthoPolyQ)::Float64 = std(x,opq.op)
std(x::Vector{Float64},mop::MultiOrthoPoly)::Float64 = sqrt(var(x,mop))

function moment(n::Int64,x::Vector{Float64},op::OrthoPoly)
    @assert 1<=n<=2 "for moments greater than 2 a quadrature rule is required"
    n==1 ? mean(x,op) : (mean(x,op)^2+var(x,op))
end

# function moment(n::Int64,x::Vector{Float64},opq::OrthoPolyQ)
#     @assert n>=1 "invalid moment"
#     if n==1
#         return mean(x,opq)
#     elseif n==2
#         return mean(x,opq)^2+var(x,opq)
#     else
#         t = Tensor(n,opq)
#         # needs IterTools for an elegant solution
#         error("This is not yet implemented")
#     end
# end


text:
using Test, DelayEmbeddings

println("\nTesting delay embeddings...")

@testset "embedding" begin

    data = Dataset(rand(10001,3))
    s = data[:, 1]; N = length(s)

    @testset "standard" begin

    	@testset "D = $(D), τ = $(τ)" for D in [2,3], τ in [2,3]
    		R = embed(s, D, τ)
    		@test R[(1+τ):end, 1] == R[1:end-τ, 2]
    		@test size(R) == (length(s) - τ*(D-1), D)
    	end
    end

    @testset "weighted" begin
    	@testset "D = $(D), τ = $(τ)" for D in [2,3], τ in [2,3]

            w = 0.1
    		R1 = embed(s, D, τ)
    		R2 = embed(s, D, τ, w)
            for γ in 0:D-1
                @test  (w^γ) * R1[1, γ+1][1] == R2[1, γ+1]
                @test  (w^γ) * R1[5, γ+1][1] == R2[5, γ+1]
            end
    	end
    end

    @testset "multi-time" begin
        D = 3
        τ1 = [2, 4]
        τ2 = [4, 8]

        R0 = embed(s, D, 2)
        R1 = embed(s, D, τ1)
        R2 = embed(s, D, τ2)

        @test R1 == R0

        R2y = R2[:, 2]
        @test R2y == R0[5:end, 1]
        @test R2[:, 1] == R0[1:end-4, 1]
        @test size(R2) == (N-maximum(τ2), 3)

        @test_throws ArgumentError embed(s, 4, τ1)
    end

end

println("\nTesting generalized embedding...")
@testset "genembed" begin
    τs = (0, 2, -7)
    js = (1, 3, 2)
    ge = GeneralizedEmbedding(τs, js)
    s = Dataset(rand(20, 3))
    x, y, z = columns(s)
    τr = τrange(s, ge)
    @testset "univariate" begin
        w = rand(20)
        em = genembed(w, τs, js)
        @test em == genembed(w, τs)
        @test em[1:3, 3] == w[1:3]
        @test em[1:3, 1] == w[1+7:3+7]
        @test em[1:3, 2] == w[1+7+2:3+7+2]
    end
    @testset "multivariate" begin
        em = genembed(s, τs, js)
        @test em[1:3, 1] == x[1+7:3+7]
        @test em[1:3, 2] == z[1+7+2:3+7+2]
        @test em[1:3, 3] == y[1:3]
    end
    @testset "weighted" begin
        ws = (1, 0, 0.1)
        em = genembed(s, τs, js; ws)
        @test em[1:3, 1] == x[1+7:3+7]
        @test em[1:3, 2] == 0 .* z[1+7+2:3+7+2]
        @test em[1:3, 3] == 0.1 .* y[1:3]
    end
    @testset "weighted integer" begin
        ws = (1, 0, -0.1)
        x = collect(1:100)
        ge = GeneralizedEmbedding(τs, js, ws)
        em = genembed(x, τs, js; ws)
        @test em[1:3, 1] == x[1+7:3+7]
        @test em[1:3, 2] == 0 .* x[1+7+2:3+7+2]
        @test em[1:3, 3] == -0.1 .* x[1:3]
    end
end


text:
"""
$(SIGNATURES)

Runs a screen given the parameters specified in `setup` using the
library `lib` and applies the `processing_func` function to the result.
"""
function run_exp(setup::ScreenSetup, lib::Library, processing_func::Function; run_idx=-1, flatten_func::Function=flatten_setup)

    guides, guide_freqs_dist = construct_library(setup, lib)

    cells, cell_phenotypes = transfect(setup, lib, guides, guide_freqs_dist)

    bin_cells = select(setup, cells, cell_phenotypes, guides)

    freqs = counts_to_freqs(bin_cells, length(guides))
    # uniform for now for all bins
    # work around for dict comprehension syntax change in v0.4->v0.5
    # https://github.com/JuliaLang/Compat.jl/issues/231
    seq_depths = Dict{Symbol, Int}()
    for binname in keys(bin_cells)
        seq_depths[binname] = setup.seq_depth
    end
    raw_data = sequencing(seq_depths, guides, freqs)
    bc_counts, genes = differences_between_bins(raw_data)

    [processing_func(bc_counts, genes)...; flatten_func(setup, lib)...; run_idx]
end

flatten_setup(setup::ScreenSetup, ::Library) = as_array(setup)
flatten_both(setup::ScreenSetup, lib::Library) = [as_array(setup)...; as_array(lib)...]


text:
export adjnodes, adjnzval, degree, dependency, expci

""" degree vector of an adjacency graph """
degree(G::SparseMatrixCSC) = vec(full(sum(G,1)))

"""
 Given an adj. matrix of a graph, return index of nodes
 adjacent to node i
"""
adjnodes(A::SparseMatrixCSC,i) = A.rowval[nzrange(A,i)]

"""
 Given an adj. matrix of a graph, return values of nodes
 adjacent to node i
"""
adjnzval(A::SparseMatrixCSC,i) = A.nzval[nzrange(A,i)]

"""
probability that any of i's neighbors will be conserved (netal)
"""
function dependency(G::SparseMatrixCSC, deg=degree(G))
    dep = 1./deg
    map!(x -> ifelse(isinf(x), 0, x), dep, dep)
end

"""
 approx expected number of conserved interactions (netal)
"""
function expci(G::SparseMatrixCSC, dep=dependency(G))
    map(i -> sum(dep[adjnodes(G,i)]), 1:size(G,1))
end

"""
expci(..)[i,j] = I[i,j] = expected # of conserved edges if i is aligned to j
"""
function expci(G1::SparseMatrixCSC,G2::SparseMatrixCSC,
               deg1=degree(G1),
               deg2=degree(G2),
               expci1=expci(G1,dependency(G1,deg1)),
               expci2=expci(G2,dependency(G2,deg2));
               normalize=true)
    z = normalize ? max(maximum(deg1),maximum(deg2)) : 1
    min.(expci1./z, (expci2./z)')
end


text:
function foo(x::Int64, y       ) 42 end
function foo(x,        y::Int64) 43 end

println(foo(1, true))
println(foo(true, 2))

function bar(x, y) foo(x, y) end
println(bar(1, true))
println(bar(true, 2))

max_stars_repo_name: gabriel-doriath-dohler/pjulia

max_issues_repo_name: gabriel-doriath-dohler/pjulia

max_forks_repo_name: gabriel-doriath-dohler/pjulia


text:
using Distributed
using JLD2

@time @everywhere include("../BioClim/src/required.jl")
cd("/Users/kiristern/Documents/GitHub/SDM")
## Get & prepare data
@time @everywhere begin
    # Load data from CSV files
    df = CSV.read("/Users/kiristern/Documents/GitHub/SDM/data/pred_prey/Puma_concolor.csv", header=true)
    # Prepare data (select columns, arrange values)
    df = prepare_gbif_data(df)
    # Separate species
    taxa_occ = [df[df.species .== u,:] for u in unique(df.species)]
    # Define coordinates range
    lon_range = (-180.0, 180.0)
    lat_range = (60.0, -50.0)
end

# Get the worldclim data by their layer number
cd("/Users/kiristern/Documents/GitHub/BioClim/")

#prediction pred
@info "Extract and crop bioclim variables"
@time wc_vars = [clip(worldclim(i), taxa_occ[1]) for i in 1:19];
# Make the prediction for each layer
@info "Predictions for each layer"
@time predictions = [bioclim(wc_vars[i], taxa_occ[1]) for i in 1:length(wc_vars)];
# Make the final prediction by taking the minimum
@info "Minimum-consensus aggregation"
@time prediction1 = reduce(minimum, predictions);
# Get the threshold for NaN given a percentile
@info "Threshold estimation"
@time threshold = first(quantile(prediction1[taxa_occ[1]], [0.05]))
@info "5% threshold:\t$(round(threshold; digits=3))"
# Filter the predictions based on the threshold
@info "Final prediction filtering"
@time for i in eachindex(prediction1.grid)
    prediction1.grid[i] < threshold && (prediction1.grid[i] = NaN)
end

#prediction prey1
@info "Extract and crop bioclim variables"
@time wc_vars = [clip(worldclim(i), taxa_occ[2]) for i in 1:19];
# Make the prediction for each layer
@info "Predictions for each layer"
@time predictions = [bioclim(wc_vars[i], taxa_occ[2]) for i in 1:length(wc_vars)];
# Make the final prediction by taking the minimum
@info "Minimum-consensus aggregation"
@time prediction2 = reduce(minimum, predictions);
# Get the threshold for NaN given a percentile
@info "Threshold estimation"
@time threshold = first(quantile(prediction2[taxa_occ[2]], [0.05]))
@info "5% threshold:\t$(round(threshold; digits=3))"
# Filter the predictions based on the threshold
@info "Final prediction filtering"
@time for i in eachindex(prediction2.grid)
    prediction2.grid[i] < threshold && (prediction2.grid[i] = NaN)
end

@info "Extract and crop bioclim variables"
@time wc_vars = [clip(worldclim(i), taxa_occ[3]) for i in 1:19];
# Make the prediction for each layer
@info "Predictions for each layer"
@time predictions = [bioclim(wc_vars[i], taxa_occ[3]) for i in 1:length(wc_vars)];
# Make the final prediction by taking the minimum
@info "Minimum-consensus aggregation"
@time prediction3 = reduce(minimum, predictions);
# Get the threshold for NaN given a percentile
@info "Threshold estimation"
@time threshold = first(quantile(prediction3[taxa_occ[3]], [0.05]))
@info "5% threshold:\t$(round(threshold; digits=3))"
# Filter the predictions based on the threshold
@info "Final prediction filtering"
@time for i in eachindex(prediction3.grid)
    prediction3.grid[i] < threshold && (prediction3.grid[i] = NaN)
end

predictions = [prediction1, prediction2, prediction3]

# ## Get the worldclim data
# @time wc_vars = map(x -> worldclim(x)[lon_range, lat_range], 1:19);
#
# ## Make predictions for all species
# @time predictions = map(x -> species_bclim(x, wc_vars), taxa_occ);

## Export predictions
@save "data/predictions.jld2" predictions

# Test import
@load "data/predictions.jld2" predictions

max_stars_repo_name: kiristern/SDM

max_issues_repo_name: kiristern/SDM

max_forks_repo_name: kiristern/SDM


text:
import Base.isempty

export AbstractPolytope,
       vertices_list,
       isempty,
       minkowski_sum

"""
    AbstractPolytope{N} <: AbstractPolyhedron{N}

Abstract type for compact convex polytopic sets.

### Notes

Every concrete `AbstractPolytope` must define the following functions:
- `vertices_list(::AbstractPolytope{N})` -- return a list of all vertices

```jldoctest; setup = :(using LazySets: subtypes)
julia> subtypes(AbstractPolytope)
4-element Array{Any,1}:
 AbstractCentrallySymmetricPolytope
 AbstractPolygon
 HPolytope
 VPolytope
```

A polytope is a bounded polyhedron (see [`AbstractPolyhedron`](@ref)).
Polytopes are compact convex sets with either of the following equivalent
properties:
1. They are the intersection of a finite number of closed half-spaces.
2. They are the convex hull of finitely many vertices.
"""
abstract type AbstractPolytope{N} <: AbstractPolyhedron{N} end

isconvextype(::Type{<:AbstractPolytope}) = true

# =============================================
# Common AbstractPolytope functions
# =============================================

"""
    isbounded(P::AbstractPolytope)

Determine whether a polytopic set is bounded.

### Input

- `P` -- polytopic set

### Output

`true` (since a polytope must be bounded).
"""
function isbounded(::AbstractPolytope)
    return true
end

"""
    isempty(P::AbstractPolytope)

Determine whether a polytope is empty.

### Input

- `P` -- abstract polytope

### Output

`true` if the given polytope contains no vertices, and `false` otherwise.

### Algorithm

This algorithm checks whether the `vertices_list` of the given polytope is empty
or not.
"""
function isempty(P::AbstractPolytope)
    return isempty(vertices_list(P))
end

"""
    isuniversal(P::AbstractPolytope{N}, [witness]::Bool=false) where {N}

Check whether a polyhedron is universal.

### Input

- `P`       -- polyhedron
- `witness` -- (optional, default: `false`) compute a witness if activated

### Output

* If `witness` option is deactivated: `false`
* If `witness` option is activated: `(false, v)` where ``v ∉ P``

### Algorithm

A witness is produced using `isuniversal(H)` where `H` is the first linear
constraint of `P`.
"""
function isuniversal(P::AbstractPolytope{N}, witness::Bool=false) where {N}
    if witness
        constraints = constraints_list(P)
        if isempty(constraints)
            return (true, N[])  # special case for polytopes without constraints
        end
        return isuniversal(constraints[1], true)
    else
        return false
    end
end

# given a polytope P, apply the linear map P to each vertex of P
# it is assumed that the interface function `vertices_list(P)` is available
@inline function _linear_map_vrep(M::AbstractMatrix, P::AbstractPolytope)
    vertices = broadcast(v -> M * v, vertices_list(P))
    m = size(M, 1) # output dimension
    if m == 1
        # TODO: substitute with Interval(convex_hull(vertices)...): convex hull 1D
        return Interval(minimum(vertices)[1], maximum(vertices)[1])
    elseif m == 2
        return VPolygon(vertices)
    else
        return VPolytope(vertices)
    end
end

function _linear_map_hrep_helper(M::AbstractMatrix, P::AbstractPolytope,
                                 algo::AbstractLinearMapAlgorithm)
    constraints = _linear_map_hrep(M, P, algo)
    m = size(M, 1) # output dimension
    if m == 1
        # TODO: create interval directly ?
        return convert(Interval, HPolytope(constraints))
    elseif m == 2
        return HPolygon(constraints)
    else
        return HPolytope(constraints)
    end
end

# the "backend" argument is ignored, used for dispatch
function _vertices_list(P::AbstractPolytope, backend)
    return vertices_list(P)
end