text:
# In-place update of the initial state distribution.
function update_a!(a::AbstractVector, α::AbstractMatrix, β::AbstractMatrix)
    @argcheck size(α, 1) == size(β, 1)
    @argcheck size(α, 2) == size(β, 2) == size(a, 1)

    K = length(a)
    c = 0.0

    for i in OneTo(K)
        a[i] = α[1,i] * β[1,i]
        c += a[i]
    end

    for i in OneTo(K)
        a[i] /= c
    end
end

# In-place update of the transition matrix.
function update_A!(A::AbstractMatrix, ξ::AbstractArray, α::AbstractMatrix, β::AbstractMatrix, LL::AbstractMatrix)
    @argcheck size(α, 1) == size(β, 1) == size(LL, 1) == size(ξ, 1)
    @argcheck size(α, 2) == size(β, 2) == size(LL, 2) == size(A, 1) == size(A, 2) == size(ξ, 2) == size(ξ, 3)

    T, K = size(LL)

    @inbounds for t in OneTo(T - 1)
        m = vec_maximum(view(LL, t+1, :))
        c = 0.0

        for i in OneTo(K), j in OneTo(K)
            ξ[t,i,j] = α[t,i] * A[i,j] * exp(LL[t + 1,j] - m) * β[t + 1,j]
            c += ξ[t,i,j]
        end

        for i in OneTo(K), j in OneTo(K)
            ξ[t,i,j] /= c
        end
    end

    fill!(A, 0.0)

    @inbounds for i in OneTo(K)
        c = 0.0

        for j in OneTo(K)
            for t in OneTo(T - 1)
                A[i,j] += ξ[t,i,j]
            end
            c += A[i,j]
        end

        for j in OneTo(K)
            A[i,j] /= c
        end
    end
end

# In-place update of the observations distributions.
function update_B!(B::AbstractVector, γ::AbstractMatrix, observations, estimator)
    @argcheck size(γ, 1) == size(observations, 1)
    @argcheck size(γ, 2) == size(B, 1)
    K = length(B)
    for i in OneTo(K)
        if sum(γ[:,i]) > 0
            B[i] = estimator(typeof(B[i]), permutedims(observations), γ[:,i])
        end
    end
end

function fit_mle!(hmm::AbstractHMM, observations; display = :none, maxiter = 100, tol=1e-3, robust = false, estimator=fit_mle)
    @argcheck display in [:none, :iter, :final]
    @argcheck maxiter >= 0

    T, K = size(observations, 1), size(hmm, 1)
    history = EMHistory(false, 0, [])

    # Allocate memory for in-place updates
    c = zeros(T)
    α = zeros(T, K)
    β = zeros(T, K)
    γ = zeros(T, K)
    ξ = zeros(T, K, K)
    LL = zeros(T, K)

    loglikelihoods!(LL, hmm, observations)
    robust && replace!(LL, -Inf => nextfloat(-Inf), Inf => log(prevfloat(Inf)))

    forwardlog!(α, c, hmm.a, hmm.A, LL)
    backwardlog!(β, c, hmm.a, hmm.A, LL)
    posteriors!(γ, α, β)

    logtot = sum(log.(c))
    (display == :iter) && println("Iteration 0: logtot = $logtot")

    for it in 1:maxiter
        update_a!(hmm.a, α, β)
        update_A!(hmm.A, ξ, α, β, LL)
        update_B!(hmm.B, γ, observations, estimator)

        # Ensure the "connected-ness" of the states,
        # this prevents case where there is no transitions
        # between two extremely likely observations.
        robust && (hmm.A .+= eps())

        @check isprobvec(hmm.a)
        @check istransmat(hmm.A)

        loglikelihoods!(LL, hmm, observations)
        robust && replace!(LL, -Inf => nextfloat(-Inf), Inf => log(prevfloat(Inf)))

        forwardlog!(α, c, hmm.a, hmm.A, LL)
        backwardlog!(β, c, hmm.a, hmm.A, LL)
        posteriors!(γ, α, β)

        logtotp = sum(log.(c))
        (display == :iter) && println("Iteration $it: logtot = $logtotp")

        push!(history.logtots, logtotp)
        history.iterations += 1

        if abs(logtotp - logtot) < tol
            (display in [:iter, :final]) && println("EM converged in $it iterations, logtot = $logtotp")
            history.converged = true
            break
        end

        logtot = logtotp
    end

    if !history.converged
        if display in [:iter, :final]
            println("EM has not converged after $(history.iterations) iterations, logtot = $logtot")
        end
    end

    history
end

mutable struct EMHistory
    converged::Bool
    iterations::Int
    logtots::Vector{Float64}
end

max_stars_repo_name: Shu-Lea-Lai/module-MultiHMMBase


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

using Markdown
using InteractiveUtils

# ╔═╡ 6fd616f6-1ea1-11eb-3814-8bfb4a096c49
using Pkg, DrWatson

# ╔═╡ 6ffe2628-1ea1-11eb-24ea-57f985146a72
begin
	@quickactivate "StatisticsWithJuliaPlutoNotebooks"
	using Distributions, Random, StatsBase, DataFrames
	using Plots, LaTeXStrings
	using Calculus, SpecialFunctions
	Random.seed!(1)
end;

# ╔═╡ ed174bc4-1ea0-11eb-1e2f-a32874cec549
md"## Listing 3.27"

# ╔═╡ 70179b62-1ea1-11eb-13ce-7fb4d54bad4b
begin
	xGrid = -5:0.01:5

	PhiA(x) = 0.5*(1+erf(x/sqrt(2)))
	PhiB(x) = cdf(Normal(),x)

	Text("Maximum difference between two CDF implementations: $(maximum(PhiA.(xGrid) - PhiB.(xGrid)))")
end

# ╔═╡ 28860946-21c1-11eb-228d-09798b3886ad
begin
	normalDensity(z) = pdf(Normal(),z)

	d0 = normalDensity.(xGrid)
	d1 = derivative.(normalDensity,xGrid)
	d2 = second_derivative.(normalDensity, xGrid)

	plot(xGrid, [d0 d1 d2], c=[:blue :red :green],label=[L"f(x)" L"f'(x)" L"f''(x)"])
	plot!([-5,5],[0,0],  color=:black, lw=0.5, xlabel="x", xlims=(-5,5), label="")
end

# ╔═╡ 70182cd8-1ea1-11eb-094d-8d0d49cf15f3
md"## End of listing 3.27"

# ╔═╡ Cell order:
# ╟─ed174bc4-1ea0-11eb-1e2f-a32874cec549
# ╠═6fd616f6-1ea1-11eb-3814-8bfb4a096c49
# ╠═6ffe2628-1ea1-11eb-24ea-57f985146a72
# ╠═70179b62-1ea1-11eb-13ce-7fb4d54bad4b
# ╠═28860946-21c1-11eb-228d-09798b3886ad
# ╟─70182cd8-1ea1-11eb-094d-8d0d49cf15f3

ext: jl


text:
#-----------------------------------------------------------------------------------------------------
# More construction method for Operation
#-----------------------------------------------------------------------------------------------------
export onsite_operation
"""
    onsite_operation(mats)

Return operation of the product of onsite operators on each site.
"""
function onsite_operation(mats::Vector{<:Matrix})
    base = size(mats[1], 1)
    len = length(mats)
    inds = [[i] for i in 1:len]
    operation(mats, inds, len, base=base)
end
#-----------------------------------------------------------------------------------------------------
"""
    onsite_operation(mats, ind, len)

Return operation of the product of onsite operator on a set of given sites.
The sites that is acted on is given by the input `ind`.
`len` indicate the size of the system.
"""
function onsite_operation(
    mats::Vector{<:Matrix},
    ind::AbstractVector{<:Integer},
    len::Int64
)
    base = size(mats[1], 1)
    inds = [[i] for i in ind]
    operation(mats, inds, len, base=base)
end

#-----------------------------------------------------------------------------------------------------
"""
    onsite_operation(mat, len)

Return operation of the product of onsite operator `mat` on each site.
`len` indicate the size of the system.
"""
function onsite_operation(mat::Matrix, len::Int64)
    mats = fill(mat, len)
    base = size(mat, 1)
    inds = [[i] for i in 1:len]
    operation(mats, inds, len, base=base)
end

#-----------------------------------------------------------------------------------------------------
export trans_inv_operation
"""
    trans_inv_operation(mat, ind, len)

Return a 1d translational invariant operator.
`mat` is the operator acting on sites `ind`, and all `ind`+i sites.
`len` indicate the system length.
"""
function trans_inv_operation(
    mat::Matrix,
    ind::AbstractVector{<:Integer},
    len::Int64
)
    mats = fill(mat, len)
    inds = [mod.(ind .+ (i-1), len) .+ 1 for i=0:len-1]
    operation(mats, inds, len)
end


text:
using LazyArrays, FillArrays, LinearAlgebra, StaticArrays, Test
import LazyArrays: MemoryLayout, DenseColumnMajor, PaddedLayout, materialize!,
                    MulAdd, Applied, ApplyLayout, arguments, DefaultApplyStyle, sub_materialize

@testset "concat" begin
    @testset "Vcat" begin
        A = @inferred(Vcat(Vector(1:10), Vector(1:20)))
        @test eltype(A) == Int
        @test @inferred(axes(A)) == (Base.OneTo(30),)
        @test @inferred(A[5]) == A[15] == 5
        @test_throws BoundsError A[31]
        @test reverse(A) == Vcat(Vector(reverse(1:20)), Vector(reverse(1:10)))
        b = Array{Int}(undef, 31)
        @test_throws DimensionMismatch copyto!(b, A)
        b = Array{Int}(undef, 30)
        @test @allocated(copyto!(b, A)) == 0
        @test b == vcat(A.args...)
        @test copy(A) isa Vcat
        @test copy(A) == A
        @test copy(A) !== A
        @test vec(A) === A
        @test A' == transpose(A) == Vector(A)'
        @test permutedims(A) == permutedims(Vector(A))

        A = @inferred(Vcat(1:10, 1:20))
        @test @inferred(length(A)) == 30
        @test @inferred(A[5]) == A[15] == 5
        @test_throws BoundsError A[31]
        @test reverse(A) == Vcat(reverse(1:20), reverse(1:10))
        b = Array{Int}(undef, 31)
        @test_throws DimensionMismatch copyto!(b, A)
        b = Array{Int}(undef, 30)
        @test_broken @allocated(copyto!(b, A)) == 0
        @test @allocated(copyto!(b, A)) ≤ 200
        @test b == vcat(A.args...)
        @test copy(A) === A
        @test vec(A) === A
        @test A' == transpose(A) == Vector(A)'
        @test A' === Hcat((1:10)', (1:20)')
        @test transpose(A) === Hcat(transpose(1:10), transpose(1:20))
        @test permutedims(A) == permutedims(Vector(A))

        A = Vcat(randn(2,10), randn(4,10))
        @test @inferred(length(A)) == 60
        @test @inferred(size(A)) == (6,10)
        @test_throws BoundsError A[61]
        @test_throws BoundsError A[7,1]
        b = Array{Float64}(undef, 7,10)
        @test_throws DimensionMismatch copyto!(b, A)
        b = Array{Float64}(undef, 6,10)
        @test_broken @allocated(copyto!(b, A)) == 0
        VERSION ≥ v"1.2" && @test @allocated(copyto!(b, A)) ≤ 200
        @test b == vcat(A.args...)
        @test copy(A) isa Vcat
        @test copy(A) == A
        @test copy(A) !== A
        @test vec(A) == vec(Matrix(A))
        @test A' == transpose(A) == Matrix(A)'
        @test permutedims(A) == permutedims(Matrix(A))

        A = Vcat(randn(2,10).+im.*randn(2,10), randn(4,10).+im.*randn(4,10))
        @test eltype(A) == ComplexF64
        @test @inferred(length(A)) == 60
        @test @inferred(size(A)) == (6,10)
        @test_throws BoundsError A[61]
        @test_throws BoundsError A[7,1]
        b = Array{ComplexF64}(undef, 7,10)
        @test_throws DimensionMismatch copyto!(b, A)
        b = Array{ComplexF64}(undef, 6,10)
        @test_broken @allocated(copyto!(b, A)) == 0
        VERSION ≥ v"1.2" && @test @allocated(copyto!(b, A)) ≤ 200
        @test b == vcat(A.args...)
        @test copy(A) isa Vcat
        @test copy(A) == A
        @test copy(A) !== A
        @test vec(A) == vec(Matrix(A))
        @test A' == Matrix(A)'
        @test transpose(A) == transpose(Matrix(A))
        @test permutedims(A) == permutedims(Matrix(A))

        @test Vcat() isa Vcat{Any,1,Tuple{}}

        A = Vcat(1,zeros(3,1))
        @test A isa AbstractMatrix
        @test A[1,1] == 1.0
        @test A[2,1] == 0.0
        @test axes(A) == (Base.OneTo(4),Base.OneTo(1))
        @test permutedims(A) == permutedims(Matrix(A))
    end
    @testset "Hcat" begin
        A = @inferred(Hcat(1:10, 2:11))
        @test_throws BoundsError A[1,3]
        @test @inferred(size(A)) == (10,2)
        @test @inferred(A[5]) == @inferred(A[5,1]) == 5
        @test @inferred(A[11]) == @inferred(A[1,2]) == 2
        b = Array{Int}(undef, 11, 2)
        @test_throws DimensionMismatch copyto!(b, A)
        b = Array{Int}(undef, 10, 2)
        @test_broken @allocated(copyto!(b, A)) == 0
        @test @allocated(copyto!(b, A)) ≤ 200
        @test b == hcat(A.args...)
        @test copy(A) === A
        @test vec(A) == vec(Matrix(A))
        @test vec(A) === Vcat(1:10,2:11)
        @test A' == Matrix(A)'
        @test A' === Vcat((1:10)', (2:11)')

        A = Hcat(Vector(1:10), Vector(2:11))
        b = Array{Int}(undef, 10, 2)
        copyto!(b, A)
        @test b == hcat(A.args...)
        @test_broken @allocated(copyto!(b, A)) == 0
        @test @allocated(copyto!(b, A)) ≤ 100
        @test copy(A) isa Hcat
        @test copy(A) == A
        @test copy(A) !== A
        @test vec(A) == vec(Matrix(A))
        @test vec(A) === Vcat(A.args...)
        @test A' == Matrix(A)'

        A = @inferred(Hcat(1, zeros(1,5)))
        @test A == hcat(1, zeros(1,5))
        @test vec(A) == vec(Matrix(A))
        @test A' == Matrix(A)'

        A = @inferred(Hcat(Vector(1:10), randn(10, 2)))
        b = Array{Float64}(undef, 10, 3)
        copyto!(b, A)
        @test b == hcat(A.args...)
        VERSION ≥ v"1.2" && @test @allocated(copyto!(b, A)) == 0
        @test vec(A) == vec(Matrix(A))

        A = Hcat(randn(5).+im.*randn(5), randn(5,2).+im.*randn(5,2))
        b = Array{ComplexF64}(undef, 5, 3)
        copyto!(b, A)
        @test b == hcat(A.args...)
        VERSION ≥ v"1.2" && @test @allocated(copyto!(b, A)) == 0
        @test vec(A) == vec(Matrix(A))
        @test A' == Matrix(A)'
        @test transpose(A) == transpose(Matrix(A))

        @testset "getindex bug" begin
            A = randn(3,3)
            H = Hcat(A,A)
            @test H[1,1] == applied(hcat,A,A)[1,1] == A[1,1]
        end
    end

    @testset "DefaultApplyStyle" begin
        v = Applied{DefaultApplyStyle}(vcat, (1, zeros(3)))
        @test v[1] == 1
        v = Applied{DefaultApplyStyle}(vcat, (1, zeros(3,1)))
        @test v[1,1] == 1
        H = Applied{DefaultApplyStyle}(hcat, (1, zeros(1,3)))
        @test H[1,1] == 1
    end

    @testset "Special pads" begin
        A = Vcat([1,2,3], Zeros(7))
        B = Vcat([1,2], Zeros(8))

        C = @inferred(A+B)
        @test C isa Vcat{Float64,1}
        @test C.args[1] isa Vector{Float64}
        @test C.args[2] isa Zeros{Float64}
        @test C == Vector(A) + Vector(B)

        B = Vcat([1,2], Ones(8))

        C = @inferred(A+B)
        @test C isa Vcat{Float64,1}
        @test C.args[1] isa Vector{Float64}
        @test C.args[2] isa Ones{Float64}
        @test C == Vector(A) + Vector(B)

        B = Vcat([1,2], randn(8))

        C = @inferred(A+B)
        @test C isa BroadcastArray{Float64}
        @test C == Vector(A) + Vector(B)

        B = Vcat(SVector(1,2), Ones(8))
        C = @inferred(A+B)
        @test C isa Vcat{Float64,1}
        @test C.args[1] isa Vector{Float64}
        @test C.args[2] isa Ones{Float64}
        @test C == Vector(A) + Vector(B)

        A = Vcat(SVector(3,4), Zeros(8))
        B = Vcat(SVector(1,2), Ones(8))
        C = @inferred(A+B)
        @test C isa Vcat{Float64,1}
        @test C.args[1] isa SVector{2,Int}
        @test C.args[2] isa Ones{Float64}
        @test C == Vector(A) + Vector(B)
    end

    @testset "Empty Vcat" begin
        @test @inferred(Vcat{Int}([1])) == [1]
        @test @inferred(Vcat{Int}()) == Int[]
    end

    @testset "in" begin
        @test 1 in Vcat(1, 1:10_000_000_000)
        @test 100_000_000 in Vcat(1, 1:10_000_000_000)
    end

    @testset "convert" begin
        for T in (Float32, Float64, ComplexF32, ComplexF64)
            Z = Vcat(zero(T),Zeros{T}(10))
            @test convert(AbstractArray,Z) ≡ Z
            @test convert(AbstractArray{T},Z) ≡ AbstractArray{T}(Z) ≡ Z
            @test convert(AbstractVector{T},Z) ≡ AbstractVector{T}(Z) ≡ Z
        end
    end

    @testset "setindex!" begin
        x = randn(5)
        y = randn(6)
        A = Vcat(x, y, 3)
        A[1] = 1
        @test A[1] == x[1] == 1
        A[6] = 2
        @test A[6] == y[1] == 2
        @test_throws MethodError A[12] = 3
        @test_throws BoundsError A[13] = 3

        x = randn(2,2); y = randn(3,2)
        A = Vcat(x,y)
        A[1,1] = 1
        @test A[1,1] == x[1,1] == 1
        A[3,1] = 2
        @test A[3,1] == y[1,1] == 2
        A[6] = 3
        @test A[1,2] == x[1,2] == 3

        x = randn(2,2); y = randn(2,3)
        B = Hcat(x,y)
        B[1,1] = 1
        @test B[1,1] == x[1,1] == 1
        B[1,3] = 2
        @test B[1,3] == y[1,1] == 2
    end

    @testset "fill!" begin
        A = Vcat([1,2,3],[4,5,6])
        fill!(A,2)
        @test A == fill(2,6)

        A = Vcat(2,[4,5,6])
        @test fill!(A,2) == fill(2,4)
        @test_throws ArgumentError fill!(A,3)

        A = Hcat([1,2,3],[4,5,6])
        fill!(A,2)
        @test A == fill(2,3,2)
    end

    @testset "Any/All" begin
        @test all(Vcat(true, Fill(true,100_000_000)))
        @test any(Vcat(false, Fill(true,100_000_000)))
        @test all(iseven, Vcat(2, Fill(4,100_000_000)))
        @test any(iseven, Vcat(2, Fill(1,100_000_000)))
        @test_throws TypeError all(Vcat(1))
        @test_throws TypeError any(Vcat(1))
    end

    @testset "isbitsunion #45" begin
        @test copyto!(Vector{Vector{Int}}(undef,6), Vcat([[1], [2], [3]], [[1], [2], [3]])) ==
            [[1], [2], [3], [1], [2], [3]]

        a = Vcat{Union{Float64,UInt8}}([1.0], [UInt8(1)])
        @test Base.isbitsunion(eltype(a))
        r = Vector{Union{Float64,UInt8}}(undef,2)
        @test copyto!(r, a) == a
        @test r == a
        @test copyto!(Vector{Float64}(undef,2), a) == [1.0,1.0]
    end

    @testset "Mul" begin
        A = Hcat([1.0 2.0],[3.0 4.0])
        B = Vcat([1.0,2.0],[3.0,4.0])

        @test MemoryLayout(typeof(A)) isa ApplyLayout{typeof(hcat)}
        @test MemoryLayout(typeof(B)) isa ApplyLayout{typeof(vcat)}
        @test A*B == Matrix(A)*Vector(B) == mul!(Vector{Float64}(undef,1),A,B) == (Vector{Float64}(undef,1) .= @~ A*B)
        @test materialize!(MulAdd(1.1,A,B,2.2,[5.0])) == 1.1*Matrix(A)*Vector(B)+2.2*[5.0]

        A = Hcat([1.0 2.0; 3 4],[3.0 4.0; 5 6])
        B = Vcat([1.0,2.0],[3.0,4.0])
        @test MemoryLayout(typeof(A)) isa ApplyLayout{typeof(hcat)}
        @test MemoryLayout(typeof(B)) isa ApplyLayout{typeof(vcat)}
        @test A*B == Matrix(A)*Vector(B) == mul!(Vector{Float64}(undef,2),A,B) == (Vector{Float64}(undef,2) .= @~ A*B)
        @test materialize!(MulAdd(1.1,A,B,2.2,[5.0,6])) ≈ 1.1*Matrix(A)*Vector(B)+2.2*[5.0,6]

        A = Hcat([1.0 2.0; 3 4],[3.0 4.0; 5 6])
        B = Vcat([1.0 2.0; 3 4],[3.0 4.0; 5 6])
        @test MemoryLayout(typeof(A)) isa ApplyLayout{typeof(hcat)}
        @test MemoryLayout(typeof(B)) isa ApplyLayout{typeof(vcat)}
        @test A*B == Matrix(A)*Matrix(B) == mul!(Matrix{Float64}(undef,2,2),A,B) == (Matrix{Float64}(undef,2,2) .= @~ A*B)
        @test materialize!(MulAdd(1.1,A,B,2.2,[5.0 6; 7 8])) ≈ 1.1*Matrix(A)*Matrix(B)+2.2*[5.0 6; 7 8]
    end

    @testset "broadcast" begin
        x = Vcat(1:2, [1,1,1,1,1], 3)
        y = 1:8
        f = (x,y) -> cos(x*y)
        @test f.(x,y) isa Vcat
        @test @inferred(broadcast(f,x,y)) == f.(Vector(x), Vector(y))

        @test (x .+ y) isa Vcat
        @test (x .+ y).args[1] isa AbstractRange
        @test (x .+ y).args[end] isa Int

        z = Vcat(1:2, [1,1,1,1,1], 3)
        @test (x .+ z) isa BroadcastArray
        @test (x + z) isa BroadcastArray
        @test Vector( x .+ z) == Vector( x + z) == Vector(x) + Vector(z)

        # Lazy mixed with Static treats as Lazy
        s = SVector(1,2,3,4,5,6,7,8)
        @test f.(x , s) isa Vcat
        @test f.(x , s) == f.(Vector(x), Vector(s))

        # these are special cased
        @test Vcat(1, Ones(5))  + Vcat(2, Fill(2.0,5)) ≡ Vcat(3, Fill(3.0,5))
        @test Vcat(SVector(1,2,3), Ones(5))  + Vcat(SVector(4,5,6), Fill(2.0,5)) ≡
            Vcat(SVector(5,7,9), Fill(3.0,5))

        H = Hcat(1, zeros(1,10))
        @test H/2 isa Hcat
        @test 2\H isa Hcat
        @test H./Ref(2) isa Hcat
        @test Ref(2).\H isa Hcat
        @test H/2  == H./Ref(2) == 2\H == Ref(2) .\ H == [1/2 zeros(1,10)]
    end

    @testset "maximum/minimum Vcat" begin
        x = Vcat(1:2, [1,1,1,1,1], 3)
        @test maximum(x) == 3
        @test minimum(x) == 1
    end

    @testset "copyto!" begin
        a = Vcat(1, Zeros(10));
        c = cache(Zeros(11));
        @test MemoryLayout(typeof(a)) isa PaddedLayout
        @test MemoryLayout(typeof(c)) isa PaddedLayout{DenseColumnMajor}
        @test copyto!(c, a) ≡ c;
        @test c.datasize[1] == 1
        @test c == a

        a = Vcat(1:3, Zeros(10))
        c = cache(Zeros(13));
        @test MemoryLayout(typeof(a)) isa PaddedLayout
        @test MemoryLayout(typeof(c)) isa PaddedLayout{DenseColumnMajor}
        @test copyto!(c, a) ≡ c;
        @test c.datasize[1] == 3
        @test c == a

        @test dot(a,a) ≡ dot(a,c) ≡ dot(c,a) ≡ dot(c,c) ≡ 14.0

        a = Vcat(1:3, Zeros(5))
        c = cache(Zeros(13));
        @test copyto!(c, a) ≡ c;
        @test c.datasize[1] == 3
        @test c[1:8] == a

        a = cache(Zeros(13)); b = cache(Zeros(15));
        @test a ≠ b
        b = cache(Zeros(13));
        a[3] = 2; b[3] = 2; b[5]=0;
        @test a == b
    end

    @testset "norm" begin
        for a in (Vcat(1,2,Fill(5,3)), Hcat([1,2],randn(2,2)), Vcat(1,Float64[])),
            p in (-Inf, 0, 0.1, 1, 2, 3, Inf)
            @test norm(a,p) ≈ norm(Array(a),p)
        end
    end

    @testset "SubVcat" begin
        A = Vcat(1,[2,3], Fill(5,10))
        V = view(A,3:5)
        @test MemoryLayout(typeof(V)) isa ApplyLayout{typeof(vcat)}
        VERSION ≥ v"1.1" && @inferred(arguments(V))
        @test arguments(V)[1] ≡ Fill(1,0)
        @test A[parentindices(V)...] == copy(V) == Array(A)[parentindices(V)...]

        A = Vcat((1:100)', Zeros(1,100),Fill(1,2,100))
        V = view(A,:,3:5)
        @test MemoryLayout(typeof(V)) isa ApplyLayout{typeof(vcat)}
        @test A[parentindices(V)...] == copy(V) == Array(A)[parentindices(V)...]
        V = view(A,2:3,3:5)
        @test MemoryLayout(typeof(V)) isa ApplyLayout{typeof(vcat)}
        @test A[parentindices(V)...] == copy(V) == Array(A)[parentindices(V)...]

        A = Hcat(1:10, Zeros(10,10))
        V = view(A,3:5,:)
        @test MemoryLayout(typeof(V)) isa ApplyLayout{typeof(hcat)}
        @test A[parentindices(V)...] == copy(V) == Array(A)[parentindices(V)...]
        V = view(A,3:5,1:4)
        @test MemoryLayout(typeof(V)) isa ApplyLayout{typeof(hcat)}
        VERSION ≥ v"1.1" && @inferred(arguments(V))
        @test arguments(V)[1] == reshape(3:5,3,1)

        v = view(A,2,1:5)
        @test MemoryLayout(typeof(v)) isa ApplyLayout{typeof(vcat)}
        @test arguments(v) == ([2], zeros(4))
        @test A[2,1:5] == copy(v) == sub_materialize(v)
    end

    @testset "Padded subarrays" begin
        a = Vcat([1,2,3],[4,5,6])
        @test sub_materialize(view(a,2:6)) == a[2:6]
        a = Vcat([1,2,3], Zeros(10))
        c = cache(Zeros(10)); c[1:3] = 1:3;
        v = view(a,2:4)
        w = view(c,2:4);
        @test MemoryLayout(typeof(a)) isa PaddedLayout{DenseColumnMajor}
        @test MemoryLayout(typeof(v)) isa PaddedLayout{DenseColumnMajor}
        @test sub_materialize(v) == a[2:4] == sub_materialize(w)
        @test sub_materialize(v) isa Vcat
        @test sub_materialize(w) isa Vcat
        A = Vcat(Eye(2), Zeros(10,2))
        V = view(A, 1:5, 1:2)
        @test sub_materialize(V) == A[1:5,1:2]
    end
end


text:
include("./shared.jl")

## 3D tests
include("../scripts/stokes_3D/Stokes3D_ve_multixpu_perf.jl")
@reset_parallel_stencil()
indsx = Int.(ceil.(LinRange(1, length(xc), 12)))
indsy = Int.(ceil.(LinRange(1, length(yc), 12)))
indsz = Int.(ceil.(LinRange(1, length(zc), 12)))
d3d   = Dict(:X=> xc[indsx], :P=>P[indsx,indsy,indsz])

@testset "Reference-tests Stokes 3D VE perf" begin
    @test_reference "reftest-files/test_Stokes3D_ve_multixpu_perf.bson" d3d by=comp
end

@reset_parallel_stencil()


text:
"A type for matched waves."
mutable struct MatchWave{T<:AbstractFloat}
    effective_waves::Vector{EffectiveWave{T}}
    average_wave::AverageWave{T}
    x_match::Vector{T} # waves are matched between average_wave.x_match
end

"Calculates the difference between the match of MatchWave.effective_waves and MatchWave.average_wave. This can be used as a proxi for convergence. "
function match_error(m_wave::MatchWave{T}; apply_norm::Function=norm) where T<:AbstractFloat
    avg_eff = AverageWave(m_wave.x_match, m_wave.effective_waves)
    j0 = findmin(abs.(m_wave.average_wave.x .- m_wave.x_match[1]))[2]
    len = length(m_wave.x_match)

    return apply_norm(m_wave.average_wave.amplitudes[j0:end,:,:][:] - avg_eff.amplitudes[:])/len
end

function MatchWave(ω::T, medium::Medium{T}, specie::Specie{T};
        radius_multiplier::T = 1.005,
        tol::T = T(1e-5), θin::T = zero(T),
        max_size::Int = 200,
        wave_effs::Vector{EffectiveWave{T}} = [zero(EffectiveWave{T})],
        x::AbstractVector{T} = [-one(T)],
        L_match::Int = 0,
        kws...
    ) where T<:Number

    k = real(ω/medium.c)

    if maximum(abs(w.k_eff) for w in wave_effs) == zero(T)
        wave_effs = effective_waves(k, medium, [specie];
            radius_multiplier=radius_multiplier,
            extinction_rescale = false,
            tol = T(10)*tol, θin=θin,
            kws...)
    end
    hankel_order = wave_effs[1].hankel_order
    # use non-dimensional effective waves
    wave_non_effs = deepcopy(wave_effs)
    # wave_effs = deepcopy(wave_effs) # uncomment to use wave_effs need to forget its pointers
    for w in wave_non_effs
       w.k_eff = w.k_eff/k
    end

    a12k = T(2)*radius_multiplier*specie.r*k
    if first(x) == - one(T)
        # using non-dimensional wave_non_effs and a12k results in non-dimensional mesh X
        L_match, X =  x_mesh_match(wave_non_effs; a12 = a12k, tol = tol, max_size=max_size)
    else
        X = x.*k
        if L_match == 0
            L_match = Int(round(length(X)/2))
        end
    end

    avg_wave_effs = [AverageWave(X[L_match:L_match+1], w) for w in wave_non_effs]
    for i in eachindex(wave_non_effs)
        wave_non_effs[i].amplitudes = wave_non_effs[i].amplitudes / norm(avg_wave_effs[i].amplitudes[1,:,1])
    end

    J = length(collect(X)) - 1
    len = (J + 1)  * (2hankel_order + 1)

    (MM_quad,b_mat) = average_wave_system(ω, X, medium, specie; tol=tol,
        radius_multiplier=radius_multiplier, hankel_order=hankel_order, θin=θin,  kws...);
    MM_mat = reshape(MM_quad, (len, len));
    b = reshape(b_mat, (len));

    (LT_mat, ER_mat, b_eff) = match_arrays(ω, wave_non_effs, L_match, X, medium, [specie]; θin=θin, a12k=a12k);

    B = b - ER_mat*b_eff
    As = (ER_mat*LT_mat + MM_mat)\B
    As_mat = reshape(As, (J+1, 2hankel_order+1, 1))

    αs = LT_mat*As + b_eff
    # use these αs to correct the magnitude of the amplitudes of the effective waves
    # and re-dimensionalise the effective wavenumbers
    for i in eachindex(wave_effs)
        wave_non_effs[i].amplitudes = αs[i] .* wave_non_effs[i].amplitudes
        wave_non_effs[i].k_eff = k * wave_non_effs[i].k_eff
    end

    # return MatchWave(wave_effs, AverageWave(hankel_order, collect(X)./k, As_mat), collect(X[L_match:end])./k)
    return MatchWave(wave_non_effs, AverageWave(hankel_order, collect(X)./k, As_mat), collect(X[L_match:end])./k)
end

"Returns (x,L), where x[L:end] is the mesh used to match with wave_effs."
function x_mesh_match(wave_effs::Vector{EffectiveWave{T}}; kws... ) where T<:AbstractFloat
    # wave_effs[end] establishes how long X should be, while wave_effs[1] estalishes how fine the mesh should be.

   # If there is only one wave, then it doesn't make sense to extend the mesh until it decays.
   # Instead we choose, arbitrarily, a quarter of the wavelength.
   # Note, having only one wave is very unusual, but tends to happen in the very low frequency limit.
    x = (length(wave_effs) > 1) ?
        x_mesh(wave_effs[end], wave_effs[1]; kws...) :
        x_mesh(wave_effs[1]; max_x = (pi/2)/abs(cos(wave_effs[1].θ_eff)*abs(wave_effs[1].k_eff)), kws...)

    x_match = x[end]
    x_max = (length(x) < length(wave_effs)*T(1.5)) ?
         x[end] + T(1.5)*(x[2] - x[1])*length(wave_effs) :
         T(2) * x[end]

    x = 0.0:(x[2] - x[1]):x_max

    L = findmin(abs.(x .- x_match))[2]

    return L, x
end


text:
using RDExample
using Test

@testset "RDExample.jl" begin
    @test my_fun(2, 2) == 10
    @test my_fun(2, 3) == 13
    @test my_fun(1, 3) == 11
    # Write your own tests here.
end

@testset "Derivative tests" begin
    @test der_my_fun(2, 1) == 2
end


text:
## Utilities for working with the bases.

"""
    num_vectors(N::Int, K::Int)

Compute the number of vectors in a basis with `N` particles and `K` sites.
"""
num_vectors(N::Int, K::Int) = binomial(N+K-1, K-1)

# Global cache of generalized Pascal's triangles for restricted num_vectors.
const triangles = Dict{Int, Matrix{Int}}()

"""
    num_vectors(N::Int, K::Int, M::Int)

Compute the number of vectors in a basis with `N` particles and `K` sites, and
a limit of `M` particles per site.
"""
function num_vectors(N::Int, K::Int, M::Int)
    0 <= N <= M * K || return 0
    N == K == 0 && return 1

    # Create a new triangle.
    if !haskey(triangles, M)
        triangles[M] = zeros(Int, 1, 1)
        triangles[M][1, 1] = 1
    end

    # Extend an existing triangle.
    if size(triangles[M], 1) < K + 1
        t_old = triangles[M]
        K_old = size(t_old, 1) - 1
        t = zeros(Int, K+1, M*K+1)
        for k in 0:K_old
            for n in 0:(M*k)
                t[k+1, n+1] = t_old[k+1, n+1]
            end
        end
        for k in (K_old+1):K
            for n in 0:(M*k)
                for m in 0:min(M, n)
                    t[k+1, n+1] += t[k, n+1-m]
                end
            end
        end
        triangles[M] = t
    end

    triangles[M][K+1, N+1]
end

"""
    num_vectors(basis::AbstractSzbasis, N::Int, K::Int)

Compute the number of vectors in a basis with `N` particles and `K` sites, and
a limit of as many particles per site as in `basis`.
"""
num_vectors(::Szbasis, N::Int, K::Int) = num_vectors(N, K)
num_vectors(basis::RestrictedSzbasis, N::Int, K::Int) = num_vectors(N, K, basis.M)

"""
    site_max(basis::AbstractSzbasis)

Get the maximum number of particles on a site in `basis`.
"""
site_max(basis::Szbasis) = basis.N
site_max(basis::RestrictedSzbasis) = basis.M

"""
    serial_num(K::Int, N::Int, v::AbstractVector{Int})

Compute the serial number of occupation vector `v` in a basis with `K` sites
and `N` particles.
"""
function serial_num(K::Int, N::Int, v::AbstractVector{Int})
    I = 1

    for mu in 1:K
        s = 0
        for nu in (mu+1):K
            s += v[nu]
        end
        for i in 0:(v[mu]-1)
            I += num_vectors(N-s-i, mu-1)
        end
    end

    I
end

"""
    serial_num(K::Int, N::Int, M::Int, v::AbstractVector{Int})

Compute the serial number of occupation vector `v` in a basis with `K` sites
and `N` particles, and a limit of `M` particles per site.
"""
function serial_num(K::Int, N::Int, M::Int, v::AbstractVector{Int})
    I = 1

    for mu in 1:K
        s = 0
        for nu in (mu+1):K
            s += v[nu]
        end
        for i in 0:(v[mu]-1)
            I += num_vectors(N-s-i, mu-1, M)
        end
    end

    I
end

"""
    serial_num(basis::AbstractSzbasis, v::AbstractVector{Int})

Compute the serial number of occupation vector `v` in `basis`.
"""
serial_num(basis::Szbasis, v::AbstractVector{Int}) = serial_num(basis.K, basis.N, v)
serial_num(basis::RestrictedSzbasis, v::AbstractVector{Int}) = serial_num(basis.K, basis.N, basis.M, v)

"""
    serial_num(basis::AbstractSzbasis, K::Int, N::Int, v::AbstractVector{Int})

Compute the serial number of occupation vector `v` in a basis with `N`
particles and `K` sites, and a limit of as many particles per site as in
`basis`.
"""
serial_num(::Szbasis, K::Int, N::Int, v::AbstractVector{Int}) = serial_num(K, N, v)
serial_num(basis::RestrictedSzbasis, K::Int, N::Int, v::AbstractVector{Int}) = serial_num(K, N, basis.M, v)

"""
    sub_serial_num(basis::AbstractSzbasis, v::AbstractVector{Int})

Compute the serial number of the reduced occupation vector `v`, which has only
a subset of the sites present in `basis`.
"""
function sub_serial_num(::Szbasis, v::AbstractVector{Int})
    K = length(v)
    N = sum(v)

    # Only one way to have no sites.
    K >= 1 || return 1
    # Only one way to have no particles.
    N >= 1 || return 1

    # Count the zero-particle case.
    I = 1

    for n in 1:(N-1)
        I += num_vectors(n, K)
    end

    I + serial_num(K, N, v)
end

function sub_serial_num(basis::RestrictedSzbasis, v::AbstractVector{Int})
    K = length(v)
    N = sum(v)

    # Only one way to have no sites.
    K >= 1 || return 1
    # Only one way to have no particles.
    N >= 1 || return 1

    # Count the zero-particle case.
    I = 1

    for n in 1:(N-1)
        I += num_vectors(n, K, basis.M)
    end

    I + serial_num(K, N, basis.M, v)
end


text:
import Base: ==

function ==(a::GeneralBigDataEntry, b::GeneralBigDataEntry)
    return a.key == b.key &&
        a.storeKey == b.storeKey &&
        a.mimeType == b.mimeType &&
        Dates.value(a.createdTimestamp - b.createdTimestamp) < 1000 &&
        Dates.value(a.lastUpdatedTimestamp - b.lastUpdatedTimestamp) < 1000 #1 second
end

function ==(a::MongodbBigDataEntry, b::MongodbBigDataEntry)
    return a.key == b.key && a.oid == b.oid
end

function ==(a::FileBigDataEntry, b::FileBigDataEntry)
    return a.key == b.key && a.filename == b.filename
end

"""
    $(SIGNATURES)
Add Big Data Entry to a DFG variable
"""
function addBigDataEntry!(var::AbstractDFGVariable, bde::AbstractBigDataEntry)
    haskey(var.bigData,bde.key) && error("BigData entry $(bde.key) already exists in variable")
    var.bigData[bde.key] = bde
    return bde
end

"""
    $(SIGNATURES)
Add Big Data Entry to distributed factor graph.
Should be extended if DFG variable is not returned by reference.
"""
function addBigDataEntry!(dfg::AbstractDFG, label::Symbol, bde::AbstractBigDataEntry)
    return addBigDataEntry!(getVariable(dfg, label), bde)
end

addDataEntry!(x...) = addBigDataEntry!(x...)

"""
$(SIGNATURES)
Add Big Data Entry to distributed factor graph.
Should be extended if DFG variable is not returned by reference.

Example

See docs for `getDataEntryElement`.

Related

addData!, getDataEntryElement, fetchData
"""
function addDataEntry!(dfg::AbstractDFG,
                       lbl::Symbol,
                       datastore::Union{FileDataStore, InMemoryDataStore},
                       descr::Symbol,
                       mimeType::AbstractString,
                       data::Vector{UInt8} )
  #
  node = isVariable(dfg, lbl) ? getVariable(dfg, lbl) : getFactor(dfg, lbl)
  # Make a big data entry in the graph - use JSON2 to just write this
  entry = GeneralBigDataEntry(dfg, node, descr, mimeType=mimeType)
  # Set it in the store
  addBigData!(datastore, entry, data)
  # Add the entry to the graph
  addBigDataEntry!(node, entry)
end
const addData! = addDataEntry!

"""
    $SIGNATURES

Does a data entry (element) exist at `key`.
"""
hasDataEntry(var::DFGVariable, key::Symbol) = haskey(var.bigData, key)
const hasBigDataEntry = hasDataEntry

"""
    $(SIGNATURES)
Get big data entry
"""
function getBigDataEntry(var::AbstractDFGVariable, key::Symbol)
    !hasDataEntry(var, key) && (error("BigData entry $(key) does not exist in variable"); return nothing)
    return var.bigData[key]
end

function getBigDataEntry(dfg::AbstractDFG, label::Symbol, key::Symbol)
    return getBigDataEntry(getVariable(dfg, label), key)
end

"""
$SIGNATURES
Get both the entry and raw data element from datastore returning as a tuple.

Notes:
- This is the counterpart to `addDataEntry!`.
- Data is identified by the node in the DFG object `dfglabel::Symbol` as well as `datalabel::Symbol`.
- The data should have been stored along with a `entry.mimeType::String` which describes the format of the data.
- ImageMagick.jl is useful for storing images in png or jpg compressed format.

Example

```julia
# some dfg object
fg = initfg()
addVariable!(fg, :x0, IIF.ContinuousScalar) # using IncrementalInference
# FileDataStore (as example)
datastore = FileDataStore(joinLogPath(fg,"datastore"))

# now some data comes in
mydata = Dict(:soundBite => randn(Float32, 10000), :meta => "something about lazy foxes and fences.")

# store the data
addDataEntry!( fg, :x0, datastore, :SOUND_DATA, "application/json", Vector{UInt8}(JSON2.write( mydata ))  )

# get/fetch the data
entry, rawData = fetchData(fg, :x0, datastore, :SOUND_DATA)

# unpack data to original format (this must be done by the user)
@show entry.mimeType # "applicatio/json"
userFormat = JSON2.read(IOBuffer(rawData))
```

Related

addDataEntry!, addData!, fetchData, fetchDataEntryElement
"""
function getDataEntryElement(dfg::AbstractDFG,
                             dfglabel::Symbol,
                             datastore::Union{FileDataStore, InMemoryDataStore},
                             datalabel::Symbol)
  #
  vari = getVariable(dfg, dfglabel)
  if !hasDataEntry(vari, datalabel)
    @error "missing data entry $datalabel in $dfglabel"
    return nothing, nothing
  end
  entry = getBigDataEntry(vari, datalabel)
  element = getBigData(datastore, entry)
  return entry, element
end
const fetchDataEntryElement = getDataEntryElement
const fetchData = getDataEntryElement

"""
    $(SIGNATURES)
Update big data entry

DevNote
- DF, unclear if `update` verb is applicable in this case, see #404
"""
function updateBigDataEntry!(var::AbstractDFGVariable,  bde::AbstractBigDataEntry)
    !haskey(var.bigData,bde.key) && (@warn "$(bde.key) does not exist in variable, adding")
    var.bigData[bde.key] = bde
    return bde
end
function updateBigDataEntry!(dfg::AbstractDFG, label::Symbol,  bde::AbstractBigDataEntry)
    # !isVariable(dfg, label) && return nothing
    return updateBigDataEntry!(getVariable(dfg, label), bde)
end

"""
    $(SIGNATURES)
Delete big data entry from the factor graph.
Note this doesn't remove it from any data stores.

Notes:
- users responsibility to delete big data in db before deleting entry
"""
function deleteBigDataEntry!(var::AbstractDFGVariable, key::Symbol)
    bde = getBigDataEntry(var, key)
    bde == nothing && return nothing
    delete!(var.bigData, key)
    return var
end
function deleteBigDataEntry!(dfg::AbstractDFG, label::Symbol, key::Symbol)
    #users responsibility to delete big data in db before deleting entry
    !isVariable(dfg, label) && return nothing
    return deleteBigDataEntry!(getVariable(dfg, label), key)
end

function deleteBigDataEntry!(var::AbstractDFGVariable, entry::AbstractBigDataEntry)
    #users responsibility to delete big data in db before deleting entry
    return deleteBigDataEntry!(var, entry.key)
end

"""
    $(SIGNATURES)
Get big data entries, Vector{AbstractBigDataEntry}
"""
function getBigDataEntries(var::AbstractDFGVariable)
    #or should we return the iterator, Base.ValueIterator{Dict{Symbol,AbstractBigDataEntry}}?
    collect(values(var.bigData))
end
function getBigDataEntries(dfg::AbstractDFG, label::Symbol)
    !isVariable(dfg, label) && return nothing
    #or should we return the iterator, Base.ValueIterator{Dict{Symbol,AbstractBigDataEntry}}?
    getBigDataEntries(getVariable(dfg, label))
end

"""
    $(SIGNATURES)
getBigDataKeys
"""
function getBigDataKeys(var::AbstractDFGVariable)
    collect(keys(var.bigData))
end
function getBigDataKeys(dfg::AbstractDFG, label::Symbol)
    !isVariable(dfg, label) && return nothing
    getBigDataKeys(getVariable(dfg, label))
end


text:
using BenchmarkTools
using Printf
using Transducers

benchmark = true

#### Record source (iterator)

struct RecordSource
    io::IOStream
end

function Base.iterate(iter::RecordSource, state = ())
    if eof(iter.io)
        close(iter.io)
        return nothing
    end
    (parse(Int64, readline(iter.io)), ())
end

Base.IteratorSize(::Type{RecordSource}) = Base.SizeUnknown()

Base.eltype(::Type{RecordSource}) = Int64

#### Factory methods

function recordSourceFactory(filename)
    factory() = RecordSource(open(filename))
    return factory
end

####

exampleRecordSource = recordSourceFactory("src/advent/day9/example.txt")
inputRecordSource = recordSourceFactory("src/advent/day9/input.txt")

####

function part1(recordSource, preambleLength)
    function isvalid(v, arr)
        for i = 1:length(arr)-1, j = i+1:length(arr)
            arr[i] + arr[j] == v && return true
        end
        false
    end

    preamble = Array{eltype(recordSource), 1}(undef, preambleLength)
    foreach(i -> preamble[i] = iterate(recordSource)[1], 1:preambleLength)
    for value ∈ recordSource
        !isvalid(value, preamble) && return value
        deleteat!(preamble, 1)
        push!(preamble, value)
    end
    nothing
end

exampleInvalidNumber = 127
@assert part1(exampleRecordSource(), 5) == exampleInvalidNumber

@show invalidNumber = part1(inputRecordSource(), 25)
benchmark && @btime part1(inputRecordSource(), 25);

####

function part2(recordSource, invalidNumber)
    window = Array{eltype(recordSource), 1}()
    for value ∈ recordSource
        push!(window, value)
        while true
            total = sum(window)
            total == invalidNumber && return minimum(window) + maximum(window)
            total < invalidNumber && break
            deleteat!(window, 1)
        end
    end
    nothing
end

@assert part2(exampleRecordSource(), exampleInvalidNumber) == 62

@show part2(inputRecordSource(), invalidNumber)
benchmark && @btime part2(inputRecordSource(), invalidNumber);


text:
function f(x, y)
    x * y
end
println(f(1, 2))

# short-hand
g(x,y) = x * y
println(g(3, 2.1))

h(x::Int64, y::Int64) = x * y
println(h(19, 65))

# default args
i(x::Int64, y::Int64, z::Int64=10) = (x + y) * z
println(i(12, 4))
println(i(12, 4, 6))

# functions with keyword args
j(x, y; mult=10) = (x + y) * mult
println("result is $(j(2, 4))")
println("result is $(j(4, 5, mult=5.6))")

# variable argument lists
function sum(args...)
    sum = 0
    for x in args
        sum += x
    end

    sum
end
println("the sum of 1..10 is $(sum(1, 2, 3, 4, 5, 6, 7, 8, 9, 10))")

max_stars_repo_name: rpillar/JuliaExamples

max_issues_repo_name: rpillar/JuliaExamples

max_forks_repo_name: rpillar/JuliaExamples


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

"""
    cancel_order(order_id)
    cancel_order(order_id, params::Dict{String,<:Any})

 Cancels an order for an Outpost.

# Arguments
- `order_id`:  The ID of the order to cancel.

"""
function cancel_order(OrderId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "POST",
        "/orders/$(OrderId)/cancel";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function cancel_order(
    OrderId, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "POST",
        "/orders/$(OrderId)/cancel",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    create_order(line_items, outpost_identifier, payment_option)
    create_order(line_items, outpost_identifier, payment_option, params::Dict{String,<:Any})

Creates an order for an Outpost.

# Arguments
- `line_items`: The line items that make up the order.
- `outpost_identifier`:  The ID or the Amazon Resource Name (ARN) of the Outpost.
- `payment_option`: The payment option for the order.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"PaymentTerm"`: The payment terms for the order.
"""
function create_order(
    LineItems,
    OutpostIdentifier,
    PaymentOption;
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "POST",
        "/orders",
        Dict{String,Any}(
            "LineItems" => LineItems,
            "OutpostIdentifier" => OutpostIdentifier,
            "PaymentOption" => PaymentOption,
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function create_order(
    LineItems,
    OutpostIdentifier,
    PaymentOption,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "POST",
        "/orders",
        Dict{String,Any}(
            mergewith(
                _merge,
                Dict{String,Any}(
                    "LineItems" => LineItems,
                    "OutpostIdentifier" => OutpostIdentifier,
                    "PaymentOption" => PaymentOption,
                ),
                params,
            ),
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    create_outpost(name, site_id)
    create_outpost(name, site_id, params::Dict{String,<:Any})

Creates an Outpost. You can specify AvailabilityZone or AvailabilityZoneId.

# Arguments
- `name`:
- `site_id`:  The ID or the Amazon Resource Name (ARN) of the site.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"AvailabilityZone"`:
- `"AvailabilityZoneId"`:
- `"Description"`:
- `"SupportedHardwareType"`:  The type of hardware for this Outpost.
- `"Tags"`: The tags to apply to the Outpost.
"""
function create_outpost(Name, SiteId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "POST",
        "/outposts",
        Dict{String,Any}("Name" => Name, "SiteId" => SiteId);
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function create_outpost(
    Name,
    SiteId,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "POST",
        "/outposts",
        Dict{String,Any}(
            mergewith(_merge, Dict{String,Any}("Name" => Name, "SiteId" => SiteId), params)
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    create_site(name)
    create_site(name, params::Dict{String,<:Any})

 Creates a site for an Outpost.

# Arguments
- `name`:

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"Description"`:
- `"Notes"`: Additional information that you provide about site access requirements,
  electrician scheduling, personal protective equipment, or regulation of equipment materials
  that could affect your installation process.
- `"OperatingAddress"`:  The location to install and power on the hardware. This address
  might be different from the shipping address.
- `"RackPhysicalProperties"`:  Information about the physical and logistical details for
  the rack at this site. For more information about hardware requirements for racks, see
  Network readiness checklist in the Amazon Web Services Outposts User Guide.
- `"ShippingAddress"`:  The location to ship the hardware. This address might be different
  from the operating address.
- `"Tags"`:  The tags to apply to a site.
"""
function create_site(Name; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "POST",
        "/sites",
        Dict{String,Any}("Name" => Name);
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function create_site(
    Name, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "POST",
        "/sites",
        Dict{String,Any}(mergewith(_merge, Dict{String,Any}("Name" => Name), params));
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    delete_outpost(outpost_id)
    delete_outpost(outpost_id, params::Dict{String,<:Any})

Deletes the Outpost.

# Arguments
- `outpost_id`:  The ID or the Amazon Resource Name (ARN) of the Outpost.

"""
function delete_outpost(OutpostId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "DELETE",
        "/outposts/$(OutpostId)";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function delete_outpost(
    OutpostId,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "DELETE",
        "/outposts/$(OutpostId)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    delete_site(site_id)
    delete_site(site_id, params::Dict{String,<:Any})

Deletes the site.

# Arguments
- `site_id`:  The ID or the Amazon Resource Name (ARN) of the site.

"""
function delete_site(SiteId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "DELETE", "/sites/$(SiteId)"; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET
    )
end
function delete_site(
    SiteId, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "DELETE",
        "/sites/$(SiteId)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    get_catalog_item(catalog_item_id)
    get_catalog_item(catalog_item_id, params::Dict{String,<:Any})

Gets information about a catalog item.

# Arguments
- `catalog_item_id`: The ID of the catalog item.

"""
function get_catalog_item(CatalogItemId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "GET",
        "/catalog/item/$(CatalogItemId)";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function get_catalog_item(
    CatalogItemId,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "GET",
        "/catalog/item/$(CatalogItemId)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    get_order(order_id)
    get_order(order_id, params::Dict{String,<:Any})

Gets an order.

# Arguments
- `order_id`: The ID of the order.

"""
function get_order(OrderId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "GET", "/orders/$(OrderId)"; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET
    )
end
function get_order(
    OrderId, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "GET",
        "/orders/$(OrderId)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    get_outpost(outpost_id)
    get_outpost(outpost_id, params::Dict{String,<:Any})

Gets information about the specified Outpost.

# Arguments
- `outpost_id`:  The ID or the Amazon Resource Name (ARN) of the Outpost.

"""
function get_outpost(OutpostId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "GET",
        "/outposts/$(OutpostId)";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function get_outpost(
    OutpostId,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "GET",
        "/outposts/$(OutpostId)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    get_outpost_instance_types(outpost_id)
    get_outpost_instance_types(outpost_id, params::Dict{String,<:Any})

Lists the instance types for the specified Outpost.

# Arguments
- `outpost_id`:  The ID or the Amazon Resource Name (ARN) of the Outpost.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"MaxResults"`:
- `"NextToken"`:
"""
function get_outpost_instance_types(
    OutpostId; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "GET",
        "/outposts/$(OutpostId)/instanceTypes";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function get_outpost_instance_types(
    OutpostId,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "GET",
        "/outposts/$(OutpostId)/instanceTypes",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    get_site(site_id)
    get_site(site_id, params::Dict{String,<:Any})

 Gets information about the specified Outpost site.

# Arguments
- `site_id`:  The ID or the Amazon Resource Name (ARN) of the site.

"""
function get_site(SiteId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "GET", "/sites/$(SiteId)"; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET
    )
end
function get_site(
    SiteId, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "GET",
        "/sites/$(SiteId)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    get_site_address(address_type, site_id)
    get_site_address(address_type, site_id, params::Dict{String,<:Any})

 Gets the site address.

# Arguments
- `address_type`:  The type of the address you request.
- `site_id`:  The ID or the Amazon Resource Name (ARN) of the site.

"""
function get_site_address(
    AddressType, SiteId; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "GET",
        "/sites/$(SiteId)/address",
        Dict{String,Any}("AddressType" => AddressType);
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function get_site_address(
    AddressType,
    SiteId,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "GET",
        "/sites/$(SiteId)/address",
        Dict{String,Any}(
            mergewith(_merge, Dict{String,Any}("AddressType" => AddressType), params)
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    list_catalog_items()
    list_catalog_items(params::Dict{String,<:Any})

Use to create a list of every item in the catalog. Add filters to your request to return a
more specific list of results. Use filters to match an item class, storage option, or EC2
family.  If you specify multiple filters, the filters are joined with an AND, and the
request returns only results that match all of the specified filters.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"EC2FamilyFilter"`:  A filter for EC2 family options for items in the catalog.  Filter
  values are case sensitive. If you specify multiple values for a filter, the values are
  joined with an OR, and the request returns all results that match any of the specified
  values.
- `"ItemClassFilter"`:  A filter for the class of items in the catalog.  Filter values are
  case sensitive. If you specify multiple values for a filter, the values are joined with an
  OR, and the request returns all results that match any of the specified values.
- `"MaxResults"`:
- `"NextToken"`:
- `"SupportedStorageFilter"`:  A filter for the storage options of items in the catalog.
  Filter values are case sensitive. If you specify multiple values for a filter, the values
  are joined with an OR, and the request returns all results that match any of the specified
  values.
"""
function list_catalog_items(; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "GET", "/catalog/items"; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET
    )
end
function list_catalog_items(
    params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "GET",
        "/catalog/items",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    list_orders()
    list_orders(params::Dict{String,<:Any})

Create a list of the Outpost orders for your Amazon Web Services account. You can filter
your request by Outpost to return a more specific list of results.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"MaxResults"`:
- `"NextToken"`:
- `"OutpostIdentifierFilter"`:  The ID or the Amazon Resource Name (ARN) of the Outpost.
"""
function list_orders(; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "GET", "/list-orders"; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET
    )
end
function list_orders(
    params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "GET",
        "/list-orders",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    list_outposts()
    list_outposts(params::Dict{String,<:Any})

Create a list of the Outposts for your Amazon Web Services account. Add filters to your
request to return a more specific list of results. Use filters to match an Outpost
lifecycle status, Availability Zone (us-east-1a), and AZ ID (use1-az1).  If you specify
multiple filters, the filters are joined with an AND, and the request returns only results
that match all of the specified filters.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"AvailabilityZoneFilter"`:  A filter for the Availability Zone (us-east-1a) of the
  Outpost.  Filter values are case sensitive. If you specify multiple values for a filter,
  the values are joined with an OR, and the request returns all results that match any of the
  specified values.
- `"AvailabilityZoneIdFilter"`:  A filter for the AZ IDs (use1-az1) of the Outpost.  Filter
  values are case sensitive. If you specify multiple values for a filter, the values are
  joined with an OR, and the request returns all results that match any of the specified
  values.
- `"LifeCycleStatusFilter"`:  A filter for the lifecycle status of the Outpost.  Filter
  values are case sensitive. If you specify multiple values for a filter, the values are
  joined with an OR, and the request returns all results that match any of the specified
  values.
- `"MaxResults"`:
- `"NextToken"`:
"""
function list_outposts(; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "GET", "/outposts"; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET
    )
end
function list_outposts(
    params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "GET", "/outposts", params; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET
    )
end

"""
    list_sites()
    list_sites(params::Dict{String,<:Any})

Lists the sites for your Amazon Web Services account.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"MaxResults"`:
- `"NextToken"`:
"""
function list_sites(; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts("GET", "/sites"; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET)
end
function list_sites(
    params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "GET", "/sites", params; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET
    )
end

"""
    list_tags_for_resource(resource_arn)
    list_tags_for_resource(resource_arn, params::Dict{String,<:Any})

Lists the tags for the specified resource.

# Arguments
- `resource_arn`: The Amazon Resource Name (ARN) of the resource.

"""
function list_tags_for_resource(
    ResourceArn; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "GET",
        "/tags/$(ResourceArn)";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function list_tags_for_resource(
    ResourceArn,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "GET",
        "/tags/$(ResourceArn)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    tag_resource(resource_arn, tags)
    tag_resource(resource_arn, tags, params::Dict{String,<:Any})

Adds tags to the specified resource.

# Arguments
- `resource_arn`: The Amazon Resource Name (ARN) of the resource.
- `tags`: The tags to add to the resource.

"""
function tag_resource(ResourceArn, Tags; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "POST",
        "/tags/$(ResourceArn)",
        Dict{String,Any}("Tags" => Tags);
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function tag_resource(
    ResourceArn,
    Tags,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "POST",
        "/tags/$(ResourceArn)",
        Dict{String,Any}(mergewith(_merge, Dict{String,Any}("Tags" => Tags), params));
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    untag_resource(resource_arn, tag_keys)
    untag_resource(resource_arn, tag_keys, params::Dict{String,<:Any})

Removes tags from the specified resource.

# Arguments
- `resource_arn`: The Amazon Resource Name (ARN) of the resource.
- `tag_keys`: The tag keys.

"""
function untag_resource(
    ResourceArn, tagKeys; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "DELETE",
        "/tags/$(ResourceArn)",
        Dict{String,Any}("tagKeys" => tagKeys);
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function untag_resource(
    ResourceArn,
    tagKeys,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "DELETE",
        "/tags/$(ResourceArn)",
        Dict{String,Any}(mergewith(_merge, Dict{String,Any}("tagKeys" => tagKeys), params));
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    update_outpost(outpost_id)
    update_outpost(outpost_id, params::Dict{String,<:Any})

 Updates an Outpost.

# Arguments
- `outpost_id`:  The ID or the Amazon Resource Name (ARN) of the Outpost.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"Description"`:
- `"Name"`:
- `"SupportedHardwareType"`:  The type of hardware for this Outpost.
"""
function update_outpost(OutpostId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "PATCH",
        "/outposts/$(OutpostId)";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function update_outpost(
    OutpostId,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "PATCH",
        "/outposts/$(OutpostId)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    update_site(site_id)
    update_site(site_id, params::Dict{String,<:Any})

 Updates the site.

# Arguments
- `site_id`:  The ID or the Amazon Resource Name (ARN) of the site.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"Description"`:
- `"Name"`:
- `"Notes"`:  Notes about a site.
"""
function update_site(SiteId; aws_config::AbstractAWSConfig=global_aws_config())
    return outposts(
        "PATCH", "/sites/$(SiteId)"; aws_config=aws_config, feature_set=SERVICE_FEATURE_SET
    )
end
function update_site(
    SiteId, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "PATCH",
        "/sites/$(SiteId)",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    update_site_address(address, address_type, site_id)
    update_site_address(address, address_type, site_id, params::Dict{String,<:Any})

 Updates the site address.   To update a site address with an order IN_PROGRESS, you must
wait for the order to complete or cancel the order.  You can update the operating address
before you place an order at the site, or after all Outposts that belong to the site have
been deactivated.

# Arguments
- `address`:  The address for the site.
- `address_type`:  The type of the address.
- `site_id`:  The ID or the Amazon Resource Name (ARN) of the site.

"""
function update_site_address(
    Address, AddressType, SiteId; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "PUT",
        "/sites/$(SiteId)/address",
        Dict{String,Any}("Address" => Address, "AddressType" => AddressType);
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function update_site_address(
    Address,
    AddressType,
    SiteId,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return outposts(
        "PUT",
        "/sites/$(SiteId)/address",
        Dict{String,Any}(
            mergewith(
                _merge,
                Dict{String,Any}("Address" => Address, "AddressType" => AddressType),
                params,
            ),
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    update_site_rack_physical_properties(site_id)
    update_site_rack_physical_properties(site_id, params::Dict{String,<:Any})

Update the physical and logistical details for a rack at a site. For more information about
hardware requirements for racks, see Network readiness checklist in the Amazon Web Services
Outposts User Guide.  To update a rack at a site with an order of IN_PROGRESS, you must
wait for the order to complete or cancel the order.

# Arguments
- `site_id`:  The ID or the Amazon Resource Name (ARN) of the site.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"FiberOpticCableType"`:  Specify the type of fiber that you will use to attach the
  Outpost to your network.
- `"MaximumSupportedWeightLbs"`:  Specify the maximum rack weight that this site can
  support. NO_LIMIT is over 2000lbs.
- `"OpticalStandard"`: Specify the type of optical standard that you will use to attach the
  Outpost to your network. This field is dependent on uplink speed, fiber type, and distance
  to the upstream device. For more information about networking requirements for racks, see
  Network in the Amazon Web Services Outposts User Guide.     OPTIC_10GBASE_SR: 10GBASE-SR
  OPTIC_10GBASE_IR: 10GBASE-IR    OPTIC_10GBASE_LR: 10GBASE-LR    OPTIC_40GBASE_SR:
  40GBASE-SR    OPTIC_40GBASE_ESR: 40GBASE-ESR    OPTIC_40GBASE_IR4_LR4L: 40GBASE-IR (LR4L)
   OPTIC_40GBASE_LR4: 40GBASE-LR4    OPTIC_100GBASE_SR4: 100GBASE-SR4
  OPTIC_100GBASE_CWDM4: 100GBASE-CWDM4    OPTIC_100GBASE_LR4: 100GBASE-LR4
  OPTIC_100G_PSM4_MSA: 100G PSM4 MSA    OPTIC_1000BASE_LX: 1000Base-LX    OPTIC_1000BASE_SX :
  1000Base-SX
- `"PowerConnector"`:  Specify the power connector that Amazon Web Services should plan to
  provide for connections to the hardware. Note the correlation between PowerPhase and
  PowerConnector.    Single-phase AC feed    L6-30P – (common in US); 30A; single phase
  IEC309 (blue) – P+N+E, 6hr; 32 A; single phase     Three-phase AC feed    AH530P7W (red)
  – 3P+N+E, 7hr; 30A; three phase    AH532P6W (red) – 3P+N+E, 6hr; 32A; three phase
- `"PowerDrawKva"`: Specify in kVA the power draw available at the hardware placement
  position for the rack.
- `"PowerFeedDrop"`:  Specify whether the power feed comes above or below the rack.
- `"PowerPhase"`:  Specify the power option that you can provide for hardware.
  Single-phase AC feed: 200 V to 277 V, 50 Hz or 60 Hz   Three-phase AC feed: 346 V to 480 V,
  50 Hz or 60 Hz
- `"UplinkCount"`: Racks come with two Outpost network devices. Depending on the supported
  uplink speed at the site, the Outpost network devices provide a variable number of uplinks.
  Specify the number of uplinks for each Outpost network device that you intend to use to
  connect the rack to your network. Note the correlation between UplinkGbps and UplinkCount.
    1Gbps - Uplinks available: 1, 2, 4, 6, 8   10Gbps - Uplinks available: 1, 2, 4, 8, 12, 16
    40 and 100 Gbps- Uplinks available: 1, 2, 4
- `"UplinkGbps"`:  Specify the uplink speed the rack should support for the connection to
  the Region.
"""
function update_site_rack_physical_properties(
    SiteId; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "PATCH",
        "/sites/$(SiteId)/rackPhysicalProperties";
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function update_site_rack_physical_properties(
    SiteId, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()
)
    return outposts(
        "PATCH",
        "/sites/$(SiteId)/rackPhysicalProperties",
        params;
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end


text:
using DrWatson
import Images
import CSV

function load_dataset(seg_df_path::String, paper_seg_labels_path::String, watershed_path::String, dapi_path::String, prior_path::String="", pciseq_path::String="";
        name::String, min_area::Float64, min_mols_per_cell=50, paper_polygons::Bool=false, dapi::Bool=false, watershed::Bool=false, pciseq=true, prior=true,
        rescale::Bool=false, store_labels::Bool=false)
    df_spatial, gene_names = B.load_df(datadir(seg_df_path))
    res = Dict(:df => df_spatial, :gene_names => gene_names, :min_area => min_area,
        :min_mols_per_cell => min_mols_per_cell, :name => name)

    if rescale
        df_spatial.x .-= minimum(df_spatial.x);
        df_spatial.y .-= minimum(df_spatial.y);

        df_spatial.x .*= 10.0; # Specifically for Allen smFISH
        df_spatial.y .*= 10.0;
    end

    df_spatial[!, :gene_name] = gene_names[df_spatial.gene]

    paper_seg_labels = nothing
    if length(paper_seg_labels_path) > 0
        paper_seg_labels = B.load_segmentation_mask(datadir(paper_seg_labels_path));
        df_spatial[!, :cell_paper] = denserank(B.staining_value_per_transcript(df_spatial, paper_seg_labels)) .- 1;
        if store_labels
            res[:paper_labels] = paper_seg_labels
        end
    end

    # df_seg_prior = B.load_df(datadir("exp_pro/allen_smfish/baysor_prior/segmentation.csv"))[1];
    # df_spatial[!, :cell_prior] = df_seg_prior.cell;

    if watershed && length(watershed_path) > 0
        watershed_labels = B.load_segmentation_mask(datadir(watershed_path))
        df_spatial[!, :cell_watershed] = denserank(B.staining_value_per_transcript(df_spatial, watershed_labels)) .- 1
        if store_labels
            res[:watershed_labels] = watershed_labels
        end
    end

    if dapi
        res[:dapi_arr] = Float16.(Images.load(datadir(dapi_path)))
    end

    if paper_polygons && (length(paper_seg_labels_path) > 0)
        res[:paper_polys] = B.extract_polygons_from_label_grid(Matrix(paper_seg_labels[1:3:end, 1:3:end]), grid_step=3.0)
    end

    if pciseq && (length(pciseq_path) > 0)
        df_pciseq = DataFrame(CSV.File(datadir(pciseq_path)), copycols=false)
        df_spatial[!, :cell_pciseq] = df_pciseq.neighbour
    end

    if prior && (length(prior_path) > 0)
        df_prior = DataFrame(CSV.File(datadir(prior_path)), copycols=false)
        df_spatial[!, :cell_prior] = df_prior.cell
    end

    return res
end

function load_merfish(seg_df_path::String="exp_pro/merfish_moffit/baysor/segmentation.csv"; min_area::Float64=25.0, paper_polygons::Bool=false, pciseq::Bool=true, kwargs...)
    res = load_dataset(
        "exp_raw/merfish_moffit/merfish_coords_adj.csv",
        "",
        "exp_raw/merfish_moffit/dapi_merged_watershed.tif",
        "exp_raw/merfish_moffit/dapi_merged.tiff",
        "exp_pro/merfish_moffit/baysor_prior/segmentation.csv";
        paper_polygons=false, min_area=min_area, name="MERFISH", kwargs...
    )

    res[:df][!, :cell_paper] = ifelse.(ismissing.(res[:df].cell), 0, denserank(res[:df].cell));

    df_seg = B.load_df(datadir(seg_df_path))[1]
    res[:df][!, :cell] = df_seg.cell;
    res[:df][!, :confidence] = df_seg.confidence;
    res[:df][!, :cluster] = df_seg.cluster;
    res[:df][!, :ncv_color] = df_seg.ncv_color;

    if paper_polygons
        res[:paper_polys] = B.boundary_polygons(res[:df], res[:df].cell_paper, grid_step=5.0, bandwidth=5.0);
    end

    if pciseq
        df_pciseq = DataFrame!(CSV.File(datadir("exp_pro/merfish_moffit/pciseq/spots.csv")))
        res[:df] = res[:df][.!in.(res[:df].gene_name, Ref(Set(["Blank-1", "Blank-2", "Blank-3", "Blank-4", "Blank-5", "Krt90"]))),:]
        res[:df][!, :cell_pciseq] = df_pciseq.neighbour
    end

    return res
end

load_osmfish(; min_area::Float64=70.0, kwargs...) =
    load_dataset(
        "exp_pro/osmfish/baysor/segmentation.csv",
        "exp_raw/osmfish/paper_segmentation.tiff",
        "exp_raw/osmfish/nuclei_watershed.tif",
        "exp_raw/osmfish/nuclei.tif",
        "exp_pro/osmfish/baysor_prior/segmentation.csv",
        "exp_pro/osmfish/pciseq/spots.csv";
        min_area=min_area, name="osmFISH", kwargs...
    )

load_starmap1020(; min_area::Float64=25.0, kwargs...) =
    load_dataset(
        "exp_pro/starmap_vis1020/baysor/segmentation.csv",
        "exp_raw/starmap_vis1020/segmentation.tiff",
        "", "",
        "exp_pro/starmap_vis1020/baysor_prior/segmentation.csv";
        dapi=false, watershed=false,
        min_area=min_area, name="STARmap 1020", kwargs...
    )

load_allen_smfish(csv_path::String="exp_pro/allen_smfish/baysor/segmentation.csv"; min_area::Float64=25.0, kwargs...) =
    load_dataset(
        csv_path,
        "exp_raw/allen_smfish/segmentation_labels_from_json_transposed.tiff",
        "exp_raw/allen_smfish/dapi_merged_watershed.tif",
        "exp_raw/allen_smfish/dapi_merged.tiff",
        "exp_pro/allen_smfish/baysor_prior/segmentation.csv",
        "exp_pro/allen_smfish/pciseq/spots.csv";
        rescale=true, min_area=min_area, name="Allen smFISH", kwargs...
    )

function load_iss(csv_path::String="exp_pro/iss_hippo/baysor/segmentation.csv"; min_area::Float64=1.0, paper_polygons::Bool=false, kwargs...)
    res = load_dataset(
        csv_path, "",
        "exp_raw/iss_hippo/CA1DapiBoundaries_4-3_right_watershed.tif",
        "exp_raw/iss_hippo/CA1DapiBoundaries_4-3_right.tif",
        "exp_pro/iss_hippo/baysor_prior/segmentation.csv",
        "exp_pro/iss_hippo/pciseq/spots.csv";
        paper_polygons=false, min_area=min_area, min_mols_per_cell=3, name="ISS", kwargs...
    )

    res[:df][!, :cell_paper] = res[:df].parent_id;
    if paper_polygons
        res[:paper_polys] = B.boundary_polygons(res[:df], res[:df].cell_paper, grid_step=1.0, bandwidth=3.0);
    end

    return res
end

max_stars_repo_name: hiraksarkar/BaysorAnalysis

max_issues_repo_name: hiraksarkar/BaysorAnalysis

max_forks_repo_name: hiraksarkar/BaysorAnalysis


text:
"""
    ssp_spinit(nuser)

Initialise SPLINT - should be called before other
SPLINT functions.

# Arguments
- `nuser::Integer`: the number of words reserved for user storage
"""
function ssp_spinit(nuser::Integer)

    nuser = Ref{Int32}(nuser)

    @ccall ssp_spinit_(nuser::Ref{Int32})::Nothing

    nothing
end

"""
    isp_s2make(istepx, istepq)

Create a spline object in memory, and return the address.
Every istep-th grid point is taken as a node point of the
spline and the grid boundaries are always included as node
points.

# Arguments
- `istepx::Integer`: steps taken in sampling the QCDNUM x grid
- `istepq::Integer`: steps taken in sampling the QCDNUM qq grid
"""
function isp_s2make(istepx::Integer, istepq::Integer)

    istepx = Ref{Int32}(istepx)
    istepq = Ref{Int32}(istepq)

    iasp = @ccall isp_s2make_(istepx::Ref{Int32}, istepq::Ref{Int32})::Int32

    iasp[]
end

"""
    isp_sxmake(istepx)

Create a 1D x spline object in memory, and return the address.
Every istep-th grid point is taken as a node point of the
spline and the grid boundaries are always included as node
points.

# Arguments
- `istepx::Integer`: steps taken in sampling the QCDNUM x grid
"""
function isp_sxmake(istepx::Integer)

    istepx = Ref{Int32}(istepx)

    iasp = @ccall isp_sxmake_(istepx::Ref{Int32})::Int32

    iasp[]
end

"""
    isp_sqmake(istepq)

Create a 1D qq spline object in memory, and return the address.
Every istep-th grid point is taken as a node point of the
spline and the grid boundaries are always included as node
points.

# Arguments
- `istepq::Integer`: steps taken in sampling the QCDNUM qq grid
"""
function isp_sqmake(istepq::Integer)

    istepq = Ref{Int32}(istepq)

    iasp = @ccall isp_sqmake_(istepq::Ref{Int32})::Int32

    iasp[]
end

"""
    ssp_uwrite(i, val)

Write something to the reserved user space.

# Arguments
- `i::Integer`: where to write, from 1 to nuser
- `val::Float`: what to write
"""
function ssp_uwrite(i::Integer, val::Float64)

    i = Ref{Int32}(i)
    val = Ref{Float64}(val)

    @ccall ssp_uwrite_(i::Ref{Int32}, val::Ref{Float64})::Nothing

    nothing
end

"""
    dsp_uread(i)

Read something from the reserved user space.

# Arguments
- `i::Integer`: where to read, from 1 to nuser

# Returns
- `val::Float64`: what is read
"""
function dsp_uread(i::Integer)

    i = Ref{Int32}(i)

    val = @ccall dsp_uread_(i::Ref{Int32})::Float64

    val[]
end

"""
    ssp_s2fill(iasp, fun, rs)

Fill the spline object by passing a function. The function
must have the signature fun(ix::Integer, iq::Integer, first::Boolean).

# Arguments
- `iasp::Integer`: address of the spline object
- `fun::Ptr{Nothing}`: function to be splined
- `rs::Float64`: set a sqrt(s) cut - 0 for no kinematic cut
"""
function ssp_s2fill(iasp::Integer, fun::Ptr{Nothing}, rs::Float64)

    iasp = Ref{Int32}(iasp)
    rs = Ref{Float64}(rs)

    @ccall ssp_s2fill_(iasp::Ref{Int32}, fun::Ptr{Cvoid},
                       rs::Ref{Float64})::Nothing

    nothing
end

"""
    ssp_sxfill(iasp, fun, iq)

Fill the 1D x spline object by passing a function. The function
must have the signature fun(ix::Integer, iq::Integer, first::Boolean).

# Arguments
- `iasp::Integer`: address of the spline object
- `fun::Ptr{Nothing}`: function to be splined
- `iq::Integer`: fixed iq value to pass to fun
"""
function ssp_sxfill(iasp::Integer, fun::Ptr{Nothing}, iq::Integer)

    iasp = Ref{Int32}(iasp)
    iq = Ref{Int32}(iq)

    @ccall ssp_sxfill_(iasp::Ref{Int32}, fun::Ptr{Cvoid},
                       iq::Ref{Int32})::Nothing

    nothing
end

"""
    ssp_sqfill(iasp, fun, ix)

Fill the 1D qq spline object by passing a function. The function
must have the signature fun(ix::Integer, iq::Integer, first::Boolean).

# Arguments
- `iasp::Integer`: address of the spline object
- `fun::Ptr{Nothing}`: function to be splined
- `ix::Integer`: fixed ix value to pass to fun
"""
function ssp_sqfill(iasp::Integer, fun::Ptr{Nothing}, ix::Integer)

    iasp = Ref{Int32}(iasp)
    ix = Ref{Int32}(ix)

    @ccall ssp_sqfill_(iasp::Ref{Int32}, fun::Ptr{Cvoid},
                       ix::Ref{Int32})::Nothing

    nothing
end

"""
    ssp_sxf123(ia, iset, def, istf, iq)

Fast structure function input for splines over x.

# Arguments
- `ia::Integer`: address of the spline object
- `iset::Integer`: QCDNUM pdf-set index
- `def::Array{Float64}`: Array of (anti-)quark
coefficients
- `istf::Integer`: structure function index
1<=>F_L, 2<=>F_2, 3<=>xF_3, 4<=>F_L'
- `iq::Integer`: index of qq value
"""
function ssp_sxf123(ia::Integer, iset::Integer, def::Array{Float64},
                    istf::Integer, iq::Integer)

    ia = Ref{Int32}(ia)
    iset = Ref{Int32}(iset)
    istf = Ref{Int32}(istf)
    iq = Ref{Int32}(iq)

    @ccall ssp_sxf123_(ia::Ref{Int32}, iset::Ref{Int32},
                       def::Ref{Float64}, istf::Ref{Int32},
                       iq::Ref{Int32})::Nothing

    nothing
end

"""
    ssp_sq123(ia, iset, def, istf, ix)

Fast structure function input for splines over qq.

# Arguments
- `ia::Integer`: address of the spline object
- `iset::Integer`: QCDNUM pdf-set index
- `def::Array{Float64}`: Array of (anti-)quark
coefficients
- `istf::Integer`: structure function index
1<=>F_L, 2<=>F_2, 3<=>xF_3, 4<=>F_L'
- `ix::Integer`: index of x value
"""
function ssp_sqf123(ia::Integer, iset::Integer, def::Array{Float64},
                    istf::Integer, ix::Integer)

    ia = Ref{Int32}(ia)
    iset = Ref{Int32}(iset)
    istf = Ref{Int32}(istf)
    ix = Ref{Int32}(ix)

    @ccall ssp_sqf123_(ia::Ref{Int32}, iset::Ref{Int32},
                       def::Ref{Float64}, istf::Ref{Int32},
                       ix::Ref{Int32})::Nothing

    nothing
end

"""
    ssp_s2f123(ia, iset, def, istf, rs)

Fast structure function input for 2D
splines over x and qq.

# Arguments
- `ia::Integer`: address of the spline object
- `iset::Integer`: QCDNUM pdf-set index
- `def::Array{Float64}`: Array of (anti-)quark
coefficients
- `istf::Integer`: structure function index
1<=>F_L, 2<=>F_2, 3<=>xF_3, 4<=>F_L'
- `rs::Float64`: sqrt(s) cut - 0 for no kinematic cut

"""
function ssp_s2f123(ia::Integer, iset::Integer, def::Array{Float64},
                    istf::Integer, rs::Float64)

    ia = Ref{Int32}(ia)
    iset = Ref{Int32}(iset)
    istf = Ref{Int32}(istf)
    rs = Ref{Float64}(rs)

    @ccall ssp_s2f123_(ia::Ref{Int32}, iset::Ref{Int32},
                       def::Ref{Float64}, istf::Ref{Int32},
                       rs::Ref{Float64})::Nothing

    nothing
end

"""
    isp_s2user(xarr, nx, qarr, nq)

Set your own node points in the spline in case
the automatic sampling fails.

The routine will discard points outside the x-qq
evolution grid, round the remaining nodes down to
the nearest grid-point and then sort them in
ascending order, discarding equal values. Thus you
are allowed to enter un-sorted scattered arrays.

# Arguments
- `xarr::Array{Float64}`: array of x values
- `nx::Integer`: length of xarr
- `qarr::Array{Float64}`: array of qq values
- `nq::Integer`: length of qarr

# Returns
- `iasp::Integer`: address of the spline object
"""
function isp_s2user(xarr::Array{Float64}, nx::Integer,
                    qarr::Array{Float64}, nq::Integer)

    nx = Ref{Int32}(nx)
    nq = Ref{Int32}(nq)

    iasp = @ccall isp_s2user_(xarr::Ref{Float64}, nx::Ref{Int32},
                              qarr::Ref{Float64}, nq::Ref{Int32})::Int32

    iasp[]
end

"""
    isp_sxuser(xarr, nx)

Set your own node points in the 1D x spline in case
the automatic sampling fails.

# Arguments
- `xarr::Array{Float64}`: array of x values
- `nx::Integer`: length of xarr

# Returns
- `iasp::Integer`: address of the spline object
"""
function isp_sxuser(xarr::Array{Float64}, nx::Integer)

    nx = Ref{Int32}(nx)

    iasp = @ccall isp_sxuser_(xarr::Ref{Float64}, nx::Ref{Int32})::Int32

    iasp[]
end

"""
    isp_squser(qarr, nq)

Set your own node points in the 1D qq spline in case
the automatic sampling fails.

# Arguments
- `qarr::Array{Float64}`: array of qq values
- `nq::Integer`: length of qarr

# Returns
- `iasp::Integer`: address of the spline object
"""
function isp_squser(qarr::Array{Float64}, nq::Integer)

    nq = Ref{Int32}(nq)

    iasp = @ccall isp_squser_(qarr::Ref{Float64}, nq::Ref{Int32})::Int32

    iasp[]
end

"""
    isp_splinetype(ia)

Get the type of spline at address `ia`.

Possible types are:
- `-1`: x spline
- `0`: not a spline
- `1`: qq spline
- `2`: x-qq spline
"""
function isp_splinetype(ia::Integer)

    ia = Ref{Int32}(ia)

    type = @ccall isp_splinetype_(ia::Ref{Int32})::Int32

    type[]
end

"""
    ssp_splims(ia)

Get node limits of spline at address `ia`.

Here, u and v refer to the x and qq
dimensions respectively.

# Returns tuple containing
- `nu::Integer`: number of nodes in u direction
- `u1::Float64`: lower u limit
- `u2::Float64`: upper u limit
- `nv::Integer`: number of nodes in v direction
- `v1::Float64`: lower v limit
- `v2::Float64`: upper v limit
- `n::Integer`: number of active nodes below kinematic limit
"""
function ssp_splims(ia)

    ia = Ref{Int32}(ia)

    nu = Ref{Int32}()
    u1 = Ref{Float64}()
    u2 = Ref{Float64}()

    nv = Ref{Int32}()
    v1 = Ref{Float64}()
    v2 = Ref{Float64}()

    n = Ref{Int32}()

    @ccall ssp_splims_(ia::Ref{Int32}, nu::Ref{Int32}, u1::Ref{Float64},
                       u2::Ref{Float64}, nv::Ref{Int32}, v1::Ref{Float64},
                       v2::Ref{Float64}, n::Ref{Int32})::Nothing

    nu[], u1[], u2[], nv[], v1[], v2[], n[]
end

"""
    ssp_unodes(ia, n, nu)

Copy u-nodes from spline at address `ia`
to local array.

# Arguments
- `ia::Integer`: address of spline
- `n::Integer`: dimension of array to copy to
- `nu::Integer`: number of u-nodes copied

# Returns
- `array::Array{Float64}`: array of u-nodes
"""
function ssp_unodes(ia::Integer, n::Integer, nu::Integer)

    ia = Ref{Int32}(ia)
    n = Ref{Int32}(n)
    nu = Ref{Int32}(nu)

    array = Array{Float64}(undef, n[])

    @ccall ssp_unodes_(ia::Ref{Int32}, array::Ref{Float64},
                       n::Ref{Int32}, nu::Ref{Int32})::Nothing

    array
end

"""
    ssp_vnodes(ia, n, nv)

Copy v-nodes from spline at address `ia`
to local array.

# Arguments
- `ia::Integer`: address of spline
- `n::Integer`: dimension of array to copy to
- `nv::Integer`: number of v-nodes copied

# Returns
- `array::Array{Float64}`: array of v-nodes
"""
function ssp_vnodes(ia::Integer, n::Integer, nv::Integer)

    ia = Ref{Int32}(ia)
    n = Ref{Int32}(n)
    nv = Ref{Int32}(nv)

    array = Array{Float64}(undef, n[])

    @ccall ssp_vnodes_(ia::Ref{Int32}, array::Ref{Float64},
                       n::Ref{Int32}, nv::Ref{Int32})::Nothing

    array
end

"""
    ssp_nprint(ia)

Print a list of nodes and grids for
the spline at address `ia`.
"""
function ssp_nprint(ia::Integer)

    ia = Ref{Int32}(ia)

    @ccall ssp_nprint_(ia::Ref{Int32})::Nothing

    nothing
end

"""
    ssp_extrapu(ia, n)

Define the extrapolation at the kinematic
limit for a spline at address `ia`.

The extrapolation index `n` can be:
- `0`: constant
- `1`: linear
- `2`: quadratic
- `3`: cubic
"""
function ssp_extrapu(ia, n)

    ia = Ref{Int32}(ia)
    n = Ref{Int32}(n)

    @ccall ssp_extrapu_(ia::Ref{Int32}, n::Ref{Int32})::Nothing

    nothing
end

"""
    ssp_extrapv(ia, n)

Define the extrapolation at the kinematic
limit for a spline at address `ia`.

The extrapolation index `n` can be:
- `0`: constant
- `1`: linear
- `2`: quadratic
- `3`: cubic
"""
function ssp_extrapv(ia, n)

    ia = Ref{Int32}(ia)
    n = Ref{Int32}(n)

    @ccall ssp_extrapv_(ia::Ref{Int32}, n::Ref{Int32})::Nothing

    nothing
end

"""
    dsp_funs1(ia, u, ichk)

Evaluate function for 1D spline.

Possible values of `ichk`:
- `-1`: extrapolate the spline
- `0`: return 0
- `1`: throw an error message

# Arguments
- `ia::Integer`: address of spline
- `u::Float64`: x or qq
- `ichk::Integer`: defines behaviour when outside
spline range
"""
function dsp_funs1(ia::Integer, u::Float64, ichk::Integer)

    ia = Ref{Int32}(ia)
    u = Ref{Float64}(u)
    ichk = Ref{Int32}(ichk)

    val = @ccall dsp_funs1_(ia::Ref{Int32}, u::Ref{Float64},
                            ichk::Ref{Int32})::Float64

    val[]
end

"""
    dsp_funs2(ia, x, q, ichk)

Evaluate function for 2D spline.

Possible values of `ichk`:
- `-1`: extrapolate the spline
- `0`: return 0
- `1`: throw an error message

# Arguments
- `ia::Integer`: address of spline
- `x::Float64`: x value
- `q::Float64`: qq value
- `ichk::Integer`: defines behaviour when outside
spline range
"""
function dsp_funs2(ia::Integer, x::Float64, q::Float64, ichk::Integer)

    ia = Ref{Int32}(ia)
    x = Ref{Float64}(x)
    q = Ref{Float64}(q)
    ichk = Ref{Int32}(ichk)

    val = @ccall dsp_funs2_(ia::Ref{Int32}, x::Ref{Float64},
                            q::Ref{Float64}, ichk::Ref{Int32})::Float64

    val[]
end

"""
    dsp_ints1(ia, u1, u2)

Evaluate integral of spline between u1 and u2.

The integration limits must lie inside the
spline range.
"""
function dsp_ints1(ia::Integer, u1::Float64, u2::Float64)

    ia = Ref{Int32}(ia)
    u1 = Ref{Float64}(u1)
    u2 = Ref{Float64}(u2)

    val = @ccall dsp_ints1_(ia::Ref{Int32}, u1::Ref{Float64},
                            u2::Ref{Float64})::Float64

    val[]
end

"""
    dsp_ints2(ia, x1, x2, q1, q2, rs, np)

Evaluate integral of spline between x1, x2, q1 and q2.
Also takes care of rscut and integrations makes use of N-point
Gauss quadrature, as defined by the choice of np.

The integration limits must lie inside the
spline range.
"""
function dsp_ints2(ia::Integer, x1::Float64, x2::Float64,
                   q1::Float64, q2::Float64, rs::Float64, np::Integer)

    ia = Ref{Int32}(ia)
    x1 = Ref{Float64}(x1)
    x2 = Ref{Float64}(x2)
    q1 = Ref{Float64}(q1)
    q2 = Ref{Float64}(q2)
    rs = Ref{Float64}(rs)
    np = Ref{Int32}(np)

    val = @ccall dsp_ints2_(ia::Ref{Int32}, x1::Ref{Float64},
                            x2::Ref{Float64}, q1::Ref{Float64},
                            q2::Ref{Float64}, rs::Ref{Float64},
                            np::Ref{Int32})::Float64

    val[]
end

"""
    ssp_erase(ia)

Clear the memory from `ia` onwards. `ia`=0 can
also be used to erase all spline objects in
memory.
"""
function ssp_erase(ia::Integer)

    ia  = Ref{Int32}(ia)

    @ccall ssp_erase_(ia::Ref{Int32})::Nothing

    nothing
end


text:
module SpeciesFeatureUtils

using MolecularGraph

# some convenience SFD constructors mapping names of species features to MolecularGraph functions...

const mg_elements = [d["Symbol"] for d in MolecularGraph.ATOMTABLE]

export sfd_names_props
const sfd_names_props = Dict(
    "hybridization" => Dict(
        :A => GraphMol,
        :compute_f => hybridization,
        :categorical => true,
        :encodable_elements => mg_elements,
        :possible_vals => [:sp3, :sp2, :sp, :none],
    ),
    "hydrogenconnected" => Dict(
        :A => GraphMol,
        :compute_f => hydrogenconnected,
        :categorical => true,
        :encodable_elements => mg_elements,
        :possible_vals => [0, 1, 2, 3], # I'm hoping this is right...
    ),
    "isaromatic" => Dict(
        :A => GraphMol,
        :compute_f => isaromatic,
        :categorical => true,
        :encodable_elements => mg_elements,
        :possible_vals => [true, false],
    ),
    "isringatom" => Dict(
        :A => GraphMol,
        :compute_f => isringatom,
        :categorical => true,
        :encodable_elements => mg_elements,
        :possible_vals => [true, false],
    ),
    "lonepair" => Dict(
        :A => GraphMol,
        :compute_f => lonepair,
        :categorical => true,
        :encodable_elements => mg_elements,
        :possible_vals => [-1, 0, 1, 2, 3],
    ),
    "pielectron" => Dict(
        :A => GraphMol,
        :compute_f => pielectron,
        :categorical => true,
        :encodable_elements => mg_elements,
        :possible_vals => [0, 1, 2], # also not certain this is correct
    ),
)

end


text:
module Hydrogen_FEM
    include("hydrogen_fem_module.jl")
    using LinearAlgebra
    using Match
    using MKL
    using .Hydrogen_FEM_module

    const NODE_TOTAL = 5000

    function construct()
        param = Hydrogen_FEM_module.Hydrogen_FEM_param(NODE_TOTAL - 1, NODE_TOTAL, 30.0, 0.0)
        val = Hydrogen_FEM_module.Hydrogen_FEM_variables(
            Symmetric(zeros(param.ELE_TOTAL, param.ELE_TOTAL)),
            Array{Float64}(undef, param.ELE_TOTAL),
            Array{Float64, 3}(undef, param.ELE_TOTAL, 2, 2),
            Array{Float64, 3}(undef, param.ELE_TOTAL, 2, 2),
            Array{Int64, 2}(undef, param.ELE_TOTAL, 2),
            Array{Float64, 2}(undef, param.ELE_TOTAL, param.ELE_TOTAL),
            Array{Float64}(undef, param.NODE_TOTAL),
            Array{Float64}(undef, param.NODE_TOTAL),
            Symmetric(zeros(param.ELE_TOTAL, param.ELE_TOTAL)))

        return param, val
    end

    function make_wavefunction(param, val)
        # データの生成
        make_data!(param, val)

        # 要素行列の生成
        make_element_matrix!(param, val)

        # 全体行列を生成
        hg_tmp, ug_tmp = make_global_matrix(param, val)

        # 境界条件処理を行う
        boundary_conditions!(param, val, hg_tmp, ug_tmp)

        # 一般化固有値問題を解く
        eigenval, phi = eigen(val.hg, val.ug)

        # 基底状態の固有ベクトルを取り出す
        val.phi = @view(phi[:,1])

        # 固有ベクトルの要素数を増やす
        resize!(val.phi, NODE_TOTAL)

        # 固有ベクトル（波動関数）を規格化
        normalize!(val)

        return eigenval[1]
    end

    rho(param, val, r) = let
        klo = 1;
        max = param.NODE_TOTAL;
        khi = max;

        # 表の中の正しい位置を二分探索で求める
        @inbounds while khi - klo > 1
            k = (khi + klo) >> 1

            if val.node_r_glo[k] > r
                khi = k
            else
                klo = k
            end
        end

        # yvec_[i] = f(xvec_[i]), yvec_[i + 1] = f(xvec_[i + 1])の二点を通る直線を代入
        tmp = (val.phi[khi] - val.phi[klo]) / (val.node_r_glo[khi] - val.node_r_glo[klo]) * (r - val.node_r_glo[klo]) + val.phi[klo]
        return tmp * tmp;
    end

    function boundary_conditions!(param, val, hg_tmp, ug_tmp)
        @inbounds for i = 1:param.ELE_TOTAL
            for j = i - 1:i + 1
                if j != 0 && j != param.NODE_TOTAL
                    # 左辺の全体行列のN行とN列を削る
                    val.hg.data[j, i] = hg_tmp.data[j, i]

                    # 右辺の全体行列のN行とN列を削る
                    val.ug.data[j, i] = ug_tmp.data[j, i]
                end
            end
        end
    end

    get_A_matrix_element(e, le, p, q) = let
        ed = float(e - 1)
        @match p begin
            1 =>
                @match q begin
                    1 => return  0.5 * le * (ed * ed + ed + 1.0 / 3.0) - le * le * (ed / 3.0 + 1.0 / 12.0)
                    2 => return -0.5 * le * (ed * ed + ed + 1.0 / 3.0) - le * le * (ed / 6.0 + 1.0 / 12.0)
                    _ => return 0.0
                end

            2 =>
                @match q begin
                    1 => return -0.5 * le * (ed * ed + ed + 1.0 / 3.0) - le * le * (ed / 6.0 + 1.0 / 12.0)
                    2 => return 0.5 * le * (ed * ed + ed + 1.0 / 3.0) - le * le * (ed / 3.0 + 0.25)
                    _ => return 0.0
                end

            _ => return 0.0
        end
    end

    get_B_matrix_element(e, le, p, q) = let
        ed = float(e - 1)
        @match p begin
            1 =>
                @match q begin
                    1 => return le * le * le * (ed * ed / 3.0 + ed / 6.0 + 1.0 / 30.0)
                    2 => return le * le * le * (ed * ed / 6.0 + ed / 6.0 + 0.05)
                    _ => return 0.0
                end

            2 =>
                @match q begin
                    1 => return le * le * le * (ed * ed / 6.0 + ed / 6.0 + 0.05)
                    2 => return le * le * le * (ed * ed / 3.0 + ed / 2.0 + 0.2)
                    _ => return 0.0
                end

            _ => return 0.0
        end
    end

    function make_data!(param, val)
        # Global節点のx座標を定義(R_MIN～R_MAX）
        dr = (param.R_MAX - param.R_MIN) / float(param.ELE_TOTAL)
        @inbounds for i = 0:param.NODE_TOTAL - 1
            # 計算領域を等分割
            val.node_r_glo[i + 1] = param.R_MIN + float(i) * dr
        end

        @inbounds for e = 1:param.ELE_TOTAL
            val.node_num_seg[e, 1] = e
            val.node_num_seg[e, 2] = e + 1
        end

        @inbounds for e = 1:param.ELE_TOTAL
            for i = 1:2
                val.node_r_ele[e, i] = val.node_r_glo[val.node_num_seg[e, i]]
            end
        end
    end

    function make_element_matrix!(param, val)
        # 各線分要素の長さを計算
        @inbounds for e = 1:param.ELE_TOTAL
            val.length[e] = abs(val.node_r_ele[e, 2] - val.node_r_ele[e, 1])
        end

        # 要素行列の各成分を計算
        @inbounds for e = 1:param.ELE_TOTAL
            le = val.length[e]
            for i = 1:2
                for j = 1:2
                    val.mat_A_ele[e, i, j] = get_A_matrix_element(e, le, i, j)
                    val.mat_B_ele[e, i, j] = get_B_matrix_element(e, le, i, j)
                end
            end
        end
    end

    function make_global_matrix(param, val)
        hg_tmp = Symmetric(zeros(param.NODE_TOTAL, param.NODE_TOTAL))
        ug_tmp = Symmetric(zeros(param.NODE_TOTAL, param.NODE_TOTAL))

        @inbounds for e = 1:param.ELE_TOTAL
            for i = 1:2
                for j = 1:2
                    hg_tmp.data[val.node_num_seg[e, i], val.node_num_seg[e, j]] += val.mat_A_ele[e, i, j]
                    ug_tmp.data[val.node_num_seg[e, i], val.node_num_seg[e, j]] += val.mat_B_ele[e, i, j]
                end
            end
        end

        return hg_tmp, ug_tmp
    end

    function normalize!(val)
        sum = 0.0
        len = length(val.phi)
        max = len - 2

        # Simpsonの公式によって数値積分する
        @inbounds @simd for i = 1:2:max
            f0 = val.phi[i] * val.phi[i] * val.node_r_glo[i] * val.node_r_glo[i]
            f1 = val.phi[i + 1] * val.phi[i + 1] * val.node_r_glo[i + 1] * val.node_r_glo[i + 1]
            f2 = val.phi[i + 2] * val.phi[i + 2] * val.node_r_glo[i + 2] * val.node_r_glo[i + 2]

            sum += (f0 + 4.0 * f1 + f2)
        end

        val.phi .*= (-1.0 / sqrt(sum * val.length[1] / 3.0))
    end
end


text:
struct Model
    f::FormulaTerm
    dict::Dict{Symbol, Any}
end

function Base.show(io::IO, m::Model)
    println(io, m.f)
    for (k, v) in m.dict
        println(io, k, ": ", v)
    end
end

"""
Capture and parse a set of expressions to generate a Model
## Arguments
* `ex::Expr`: an expression parsed as a `Formula`
* `fe::Expr`: Fixed effect expression.   You can add an arbitrary number of high dimensional fixed effects, separated with `+`.  Interact multiple categorical variables using `&` .     Interact a categorical variable with a continuous variable using `&`.   Alternative, use `*` to add a categorical variable and its interaction with a continuous variable. Variables must be of type CategoricalArray (use `categorical` to convert a variable to a `CategoricalArray`).
* `vcov::Expr`: Vcov formula. Default to `simple`. `robust` and `cluster()` are also implemented
* `weights::Expr`: Weight variable. Corresponds to analytical weights
* `subset::Expr`: Expression of the form State .>= 30

### Returns
* `::Model`: a `Model` struct

## Detail
A typical formula is composed of one dependent variable, exogeneous variables, endogeneous variables, and instruments
```
depvar ~ exogeneousvars + (endogeneousvars ~ instrumentvars
```

### Examples
```julia
using DataFrames, RDatasets, FixedEffectModels
df = dataset("plm", "Cigar")
df.StateC =  categorical(df.State)
df.YearC =  categorical(df.Year)
reg(df, @model(Sales ~ NDI, weights = Pop))
@model(Sales ~ NDI, fe = StateC, vcov = robust)
@model(Sales ~ NDI, fe = StateC + YearC, weights = Pop, vcov = cluster(StateC))

```
"""
macro model(ex, kws...)
    f = FixedEffectModels.terms!(FixedEffectModels.sort_terms!(FixedEffectModels.parse!(ex)))
    d = Dict{Symbol, Any}()
    for kw in kws
       isa(kw, Expr) &&  kw.head== :(=) || throw("All arguments of @model, except the first one, should be keyboard arguments")
       d[kw.args[1]] = kw.args[2]
    end
    :(Model($f, $d))
end


text:
# ------------------------------------------------------------------
# Licensed under the ISC License. See LICENSE in the project root.
# ------------------------------------------------------------------

"""
    ProductPartitioner(p₁, p₂)

A method for partitioning spatial objects using the product of two
partitioners `p₁` and `p₂`.
"""
struct ProductPartitioner{P1<:AbstractFunctionPartitioner,
                          P2<:AbstractFunctionPartitioner} <: AbstractFunctionPartitioner
  p₁::P1
  p₂::P2
end

(p::ProductPartitioner)(i, j) = p.p₁(i, j) * p.p₂(i, j)

Base.:*(p₁, p₂) = ProductPartitioner(p₁, p₂)