text:
# written by Poom Chiarawongse <eight1911@gmail.com>

module util

    export gini, entropy, zero_one, q_bi_sort!, hypergeometric, check_input

    function assign(Y :: AbstractVector{T}, list :: AbstractVector{T}) where T
        dict = Dict{T, Int}()
        @simd for i in 1:length(list)
            @inbounds dict[list[i]] = i
        end

        _Y = Array{Int}(undef, length(Y))
        @simd for i in 1:length(Y)
            @inbounds _Y[i] = dict[Y[i]]
        end

        return list, _Y
    end

    function assign(Y :: AbstractVector{T}) where T
        set = Set{T}()
        for y in Y
            push!(set, y)
        end
        list = collect(set)
        return assign(Y, list)
    end

    @inline function zero_one(ns, n)
        return 1.0 - maximum(ns) / n
    end

    @inline function gini(ns, n)
        s = 0.0
        @simd for k in ns
            s += k * (n - k)
        end
        return s / (n * n)
    end

    # compute table of values i*log(i) for integers in 0 <= i <= maxvalue
    # where tables[i+1] = i * log(i)
    # (0*log(0) is set to 0 for convenience when computing entropy)
    function compute_entropy_terms(maxvalue)
        entropy_terms = zeros(Float64, maxvalue+1)
        for i in 1:maxvalue
            entropy_terms[i+1] = i * log(i)
        end
        return entropy_terms
    end

    # returns the entropy of ns/n, ns is an array of integers
    # and entropy_terms are precomputed entropy terms
    @inline function entropy(ns::AbstractVector{U}, n, entropy_terms) where {U <: Integer}
        s = 0.0
        for k in ns
            s += entropy_terms[k+1]
        end
        return log(n) - s / n
    end

    @inline function entropy(ns, n)
        s = 0.0
        @simd for k in ns
            if k > 0
                s += k * log(k)
            end
        end
        return log(n) - s / n
    end

    # adapted from the Julia Base.Sort Library
    @inline function partition!(v, w, pivot, region)
        i, j = 1, length(region)
        r_start = region.start - 1
        @inbounds while true
            while w[i] <= pivot; i += 1; end;
            while w[j]  > pivot; j -= 1; end;
            i >= j && break
            ri = r_start + i
            rj = r_start + j
            v[ri], v[rj] = v[rj], v[ri]
            w[i], w[j] = w[j], w[i]
            i += 1; j -= 1
        end
        return j
    end

    @inline function partitionv2!(v, w, pivot, region)
        i, j = 1, length(region)
        r_start = 0
        @inbounds while true
            while w[i] <= pivot; i += 1; end;
            while w[j]  > pivot; j -= 1; end;
            i >= j && break
            ri = r_start + i
            rj = r_start + j
            v[ri], v[rj] = v[rj], v[ri]
            w[i], w[j] = w[j], w[i]
            i += 1; j -= 1
        end
        return j
    end

    # adapted from the Julia Base.Sort Library
    function insert_sort!(v, w, lo, hi, offset)
        @inbounds for i = lo+1:hi
            j = i
            x = v[i]
            y = w[offset+i]
            while j > lo
                if x < v[j-1]
                    v[j] = v[j-1]
                    w[offset+j] = w[offset+j-1]
                    j -= 1
                    continue
                end
                break
            end
            v[j] = x
            w[offset+j] = y
        end
        return v
    end

    @inline function _selectpivot!(v, w, lo, hi, offset)
        @inbounds begin
            mi = (lo+hi)>>>1

            # sort the values in v[lo], v[mi], v[hi]

            if v[mi] < v[lo]
                v[mi], v[lo] = v[lo], v[mi]
                w[offset+mi], w[offset+lo] = w[offset+lo], w[offset+mi]
            end
            if v[hi] < v[mi]
                if v[hi] < v[lo]
                    v[lo], v[mi], v[hi] = v[hi], v[lo], v[mi]
                    w[offset+lo], w[offset+mi], w[offset+hi] = w[offset+hi], w[offset+lo], w[offset+mi]
                else
                    v[hi], v[mi] = v[mi], v[hi]
                    w[offset+hi], w[offset+mi] = w[offset+mi], w[offset+hi]
                end
            end

            # move v[mi] to v[lo] and use it as the pivot
            v[lo], v[mi] = v[mi], v[lo]
            w[offset+lo], w[offset+mi] = w[offset+mi], w[offset+lo]
            v_piv = v[lo]
            w_piv = w[offset+lo]
        end

        # return the pivot
        return v_piv, w_piv
    end

    # adapted from the Julia Base.Sort Library
    @inline function _bi_partition!(v, w, lo, hi, offset)
        pivot, w_piv = _selectpivot!(v, w, lo, hi, offset)
        # pivot == v[lo], v[hi] > pivot
        i, j = lo, hi
        @inbounds while true
            i += 1; j -= 1
            while v[i] < pivot; i += 1; end;
            while pivot < v[j]; j -= 1; end;
            i >= j && break
            v[i], v[j] = v[j], v[i]
            w[offset+i], w[offset+j] = w[offset+j], w[offset+i]
        end
        v[j], v[lo] = pivot, v[j]
        w[offset+j], w[offset+lo] = w_piv, w[offset+j]

        # v[j] == pivot
        # v[k] >= pivot for k > j
        # v[i] <= pivot for i < j
        return j
    end

    # adapted from the Julia Base.Sort Library
    # adapted from the Julia Base.Sort Library
    # this sorts v[lo:hi] and w[offset+lo, offset+hi]
    # simultaneously by the values in v[lo:hi]
    const SMALL_THRESHOLD  = 20
    function q_bi_sort!(v, w, lo, hi, offset)
        @inbounds while lo < hi
            hi-lo <= SMALL_THRESHOLD && return insert_sort!(v, w, lo, hi, offset)
            j = _bi_partition!(v, w, lo, hi, offset)
            if j-lo < hi-j
                # recurse on the smaller chunk
                # this is necessary to preserve O(log(n))
                # stack space in the worst case (rather than O(n))
                lo < (j-1) && q_bi_sort!(v, w, lo, j-1, offset)
                lo = j+1
            else
                j+1 < hi && q_bi_sort!(v, w, j+1, hi, offset)
                hi = j-1
            end
        end
        return v
    end

    # The code function below is a small port from numpy's library
    # library which is distributed under the 3-Clause BSD license.
    # The rest of DecisionTree.jl is released under the MIT license.

    # ported by Poom Chiarawongse <eight1911@gmail.com>

    # this is the code for efficient generation
    # of hypergeometric random variables ported from numpy.random
    function hypergeometric(good, bad, sample, rng)

        @inline function loggam(x)
            x0 = x
            n = 0
            if (x == 1.0 || x == 2.0)
                return 0.0
            elseif x <= 7.0
                n = Int(floor(7 - x))
                x0 = x + n
            end
            x2 = 1.0 / (x0*x0)
            xp = 6.2831853071795864769252867665590 # Tau
            gl0 = -1.39243221690590e+00
            gl0 = gl0 * x2 + 1.796443723688307e-01
            gl0 = gl0 * x2 - 2.955065359477124e-02
            gl0 = gl0 * x2 + 6.410256410256410e-03
            gl0 = gl0 * x2 - 1.917526917526918e-03
            gl0 = gl0 * x2 + 8.417508417508418e-04
            gl0 = gl0 * x2 - 5.952380952380952e-04
            gl0 = gl0 * x2 + 7.936507936507937e-04
            gl0 = gl0 * x2 - 2.777777777777778e-03
            gl0 = gl0 * x2 + 8.333333333333333e-02
            gl = gl0/x0 + 0.5*log(xp) + (x0-0.5)*log(x0) - x0
            if x <= 7.0
                @simd for k in 1:n
                    gl -= log(x0 - k)
                end
            end
            return gl
        end

        @inline function hypergeometric_hyp(good, bad, sample)
            d1 = bad + good - sample
            d2 = min(bad, good)

            Y = d2
            K = sample
            while Y > 0
                Y -= floor(UInt, rand(rng) + Y/(d1 + K))
                K -= 1
                if K == 0
                    break
                end
            end

            Z = d2 - Y
            return if good > bad
                sample - Z
            else
                Z
            end
        end

        @inline function hypergeometric_hrua(good, bad, sample)
            mingoodbad = min(good, bad)
            maxgoodbad = max(good, bad)
            popsize = good + bad
            m = min(sample, popsize - sample)
            d4 = mingoodbad / popsize
            d5 = 1.0 - d4
            d6 = m*d4 + 0.5
            d7 = sqrt((popsize - m) * sample * d4 * d5 / (popsize - 1) + 0.5)
            # d8 = 2*sqrt(2/e) * d7 + (3 - 2*sqrt(3/e))
            d8 = 1.7155277699214135*d7 + 0.8989161620588988
            d9 = floor(UInt, (m + 1) * (mingoodbad + 1) / (popsize + 2))
            d10 = (loggam(d9+1) + loggam(mingoodbad-d9+1) + loggam(m-d9+1) +
                   loggam(maxgoodbad-m+d9+1))
            d11 = min(m+1, mingoodbad+1, floor(UInt, d6+16*d7))

            while true
                X = rand(rng)
                Y = rand(rng)
                W = d6 + d8*(Y - 0.5)/X

                (W < 0.0 || W >= d11) && continue
                Z = floor(Int, W)
                T = d10 - (loggam(Z+1) + loggam(mingoodbad-Z+1) + loggam(m-Z+1) +
                           loggam(maxgoodbad-m+Z+1))
                (X*(4.0-X)-3.0) <= T && break
                (X*(X-T) >= 1) && continue
                (2.0*log(X) <= T) && break
            end

            if good > bad
                Z = m - Z
            end

            return if m < sample
                good - Z
            else
                Z
            end
        end

        return if sample > 10
            hypergeometric_hrua(good, bad, sample)
        else
            hypergeometric_hyp(good, bad, sample)
        end
    end

    function check_input(
            X                   :: AbstractMatrix{S},
            Y                   :: AbstractVector{T},
            W                   :: AbstractVector{U},
            max_features        :: Int,
            max_depth           :: Int,
            min_samples_leaf    :: Int,
            min_samples_split   :: Int,
            min_purity_increase :: Float64) where {S, T, U}
        n_samples, n_features = size(X)
        if length(Y) != n_samples
            throw("dimension mismatch between X and Y ($(size(X)) vs $(size(Y))")
        elseif length(W) != n_samples
            throw("dimension mismatch between X and W ($(size(X)) vs $(size(W))")
        elseif max_depth < -1
            throw("unexpected value for max_depth: $(max_depth) (expected:"
                * " max_depth >= 0, or max_depth = -1 for infinite depth)")
        elseif n_features < max_features
            throw("number of features $(n_features) is less than the number "
                * "of max features $(max_features)")
        elseif max_features < 0
            throw("number of features $(max_features) must be >= zero ")
        elseif min_samples_leaf < 1
            throw("min_samples_leaf must be a positive integer "
                * "(given $(min_samples_leaf))")
        elseif min_samples_split < 2
            throw("min_samples_split must be at least 2 "
                * "(given $(min_samples_split))")
        end
    end

    function check_input(
            X                   :: AbstractMatrix{S},
            Y                   :: AbstractVector{Float64},
            W                   :: AbstractVector{Int},
            indX                :: AbstractVector{Int},
            const_mtry          :: Bool,
            mpois               :: Int,
            max_depth           :: Int,
            min_samples_leaf    :: Int,
            min_samples_split   :: Int) where {S}
        n_samples, n_features = size(X)
        if length(Y) != n_samples
            throw("dimension mismatch between X and Y ($(size(X)) vs $(size(Y))")
        elseif length(W) != n_samples
            throw("dimension mismatch between X and T ($(size(X)) vs $(size(W))")
        elseif max_depth < -1
            throw("unexpected value for max_depth: $(max_depth) (expected:"
                * " max_depth >= 0, or max_depth = -1 for infinite depth)")
        elseif mpois < -1
            throw("unexpected value for max_depth: $(mpois) (expected:"
                * " mpois >= 0, or mpois = -1 for by default value)")
        elseif min_samples_leaf < 1
            throw("min_samples_leaf must be a positive integer "
                * "(given $(min_samples_leaf))")
        elseif min_samples_split < 2
            throw("min_samples_split must be at least 2 "
                * "(given $(min_samples_split))")
        end
    end

end


text:
using JuliaEuler
using Documenter

DocMeta.setdocmeta!(JuliaEuler, :DocTestSetup, :(using JuliaEuler); recursive=true)

makedocs(;
    modules=[JuliaEuler],
    authors="Brynjar Smári Bjarnason <binni@binnisb.com> and contributors",
    repo="https://github.com/binnisb/JuliaEuler.jl/blob/{commit}{path}#{line}",
    sitename="JuliaEuler.jl",
    format=Documenter.HTML(;
        prettyurls=get(ENV, "CI", "false") == "true",
        canonical="https://binnisb.github.io/JuliaEuler.jl",
        assets=String[],
    ),
    pages=[
        "Home" => "index.md",
    ],
)

deploydocs(;
    repo="github.com/binnisb/JuliaEuler.jl",
)

ext: jl

max_stars_repo_name: binnisb/JuliaEuler

max_issues_repo_name: binnisb/JuliaEuler

max_forks_repo_name: binnisb/JuliaEuler


text:
using Test
using Dash: Resource, ResourcePkg, ResourcesRegistry,
    AppResource, AppRelativeResource, AppExternalResource, AppCustomResource, AppAssetResource, ApplicationResources,
    isdynamic, register_package!

@testset "lazy + dynamic" begin
    test_resource = Resource(
        relative_package_path = "test.js",
        dynamic = true
    )
    @test test_resource.async == :lazy
    @test isdynamic(test_resource, true)
    @test isdynamic(test_resource, false)

    test_resource = Resource(
        relative_package_path = "test.js",
        dynamic = false
    )
    @test test_resource.async == :none
    @test !isdynamic(test_resource, true)
    @test !isdynamic(test_resource, false)

    @test_throws ArgumentError test_resource = Resource(
        relative_package_path = "test.js",
        dynamic = false,
        async = :false
    )

    test_resource = Resource(
        relative_package_path = "test.js",
        async = :lazy
    )

    @test test_resource.async == :lazy
    @test isdynamic(test_resource, true)
    @test isdynamic(test_resource, false)

    test_resource = Resource(
        relative_package_path = "test.js",
        async = :eager
    )

    @test test_resource.async == :eager
    @test !isdynamic(test_resource, true)
    @test isdynamic(test_resource, false)
end

@testset "application resources base registry" begin
    test_registry = ResourcesRegistry(
        dash_dependency = (
            dev = ResourcePkg(
                "dash_renderer",
                "path",
                [
                    Resource(
                        relative_package_path = ["path1.dev.js", "path2.dev.js"],
                        external_url = ["external_path1.dev.js", "external_path2.dev.js"]
                    )
                ]
            ),
            prod = ResourcePkg(
                "dash_renderer",
                "path",
                [
                    Resource(
                        relative_package_path = ["path1.prod.js", "path2.prod.js"],
                        external_url = ["external_path1.prod.js", "external_path2.prod.js"]
                    )
                ]
            )
        ),
        dash_renderer = ResourcePkg(
            "dash_renderer",
            "path",
            [
                Resource(
                    relative_package_path = "dash-renderer/dash_renderer.min.js",
                    dev_package_path = "dash-renderer/dash_renderer.dev.js",
                    external_url = "https://dash_renderer.min.js"
                )
            ]
        )
    )

    test_app = dash(include_assets_files = false)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 0
    @test length(app_resources.js) == 3
    @test app_resources.js[1] isa AppRelativeResource
    @test app_resources.js[1].namespace == "dash_renderer"
    @test app_resources.js[1].relative_path == "path1.prod.js"

    @test app_resources.js[2] isa AppRelativeResource
    @test app_resources.js[2].namespace == "dash_renderer"
    @test app_resources.js[2].relative_path == "path2.prod.js"

    @test app_resources.js[end] isa AppRelativeResource
    @test app_resources.js[end].namespace == "dash_renderer"
    @test app_resources.js[end].relative_path == "dash-renderer/dash_renderer.min.js"

    @test haskey(app_resources.files, "dash_renderer")
    @test app_resources.files["dash_renderer"].base_path == "path"
    @test "path1.prod.js" in app_resources.files["dash_renderer"].files
    @test "path2.prod.js" in app_resources.files["dash_renderer"].files
    @test "dash-renderer/dash_renderer.min.js" in app_resources.files["dash_renderer"].files

    test_app = dash(include_assets_files = false)
    enable_dev_tools!(test_app, debug = false, dev_tools_serve_dev_bundles = true)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 0
    @test length(app_resources.js) == 3

    @test app_resources.js[1] isa AppRelativeResource
    @test app_resources.js[1].namespace == "dash_renderer"
    @test app_resources.js[1].relative_path == "path1.prod.js"

    @test app_resources.js[2] isa AppRelativeResource
    @test app_resources.js[2].namespace == "dash_renderer"
    @test app_resources.js[2].relative_path == "path2.prod.js"

    @test app_resources.js[end] isa AppRelativeResource
    @test app_resources.js[end].namespace == "dash_renderer"
    @test app_resources.js[end].relative_path == "dash-renderer/dash_renderer.dev.js"

    @test haskey(app_resources.files, "dash_renderer")
    @test app_resources.files["dash_renderer"].base_path == "path"
    @test "path1.prod.js" in app_resources.files["dash_renderer"].files
    @test "path2.prod.js" in app_resources.files["dash_renderer"].files
    @test "dash-renderer/dash_renderer.dev.js" in app_resources.files["dash_renderer"].files

    test_app = dash(include_assets_files = false)
    enable_dev_tools!(test_app, debug = false, dev_tools_props_check = true, dev_tools_serve_dev_bundles = true)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 0
    @test length(app_resources.js) == 3

    @test app_resources.js[1] isa AppRelativeResource
    @test app_resources.js[1].namespace == "dash_renderer"
    @test app_resources.js[1].relative_path == "path1.dev.js"

    @test app_resources.js[2] isa AppRelativeResource
    @test app_resources.js[2].namespace == "dash_renderer"
    @test app_resources.js[2].relative_path == "path2.dev.js"

    @test app_resources.js[end] isa AppRelativeResource
    @test app_resources.js[end].namespace == "dash_renderer"
    @test app_resources.js[end].relative_path == "dash-renderer/dash_renderer.dev.js"

    @test haskey(app_resources.files, "dash_renderer")
    @test app_resources.files["dash_renderer"].base_path == "path"
    @test "path1.dev.js" in app_resources.files["dash_renderer"].files
    @test "path2.dev.js" in app_resources.files["dash_renderer"].files
    @test "dash-renderer/dash_renderer.dev.js" in app_resources.files["dash_renderer"].files

    test_app = dash(serve_locally = false, include_assets_files = false)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 0
    @test length(app_resources.js) == 3

    @test app_resources.js[1] isa AppExternalResource
    @test app_resources.js[1].url == "external_path1.prod.js"

    @test app_resources.js[2] isa AppExternalResource
    @test app_resources.js[2].url == "external_path2.prod.js"

    @test app_resources.js[end] isa AppExternalResource
    @test app_resources.js[end].url == "https://dash_renderer.min.js"

    @test isempty(app_resources.files)
end

@testset "application resources dynamic" begin
    test_registry = ResourcesRegistry(
        dash_dependency = (
            dev = ResourcePkg(
                "dash_renderer",
                "path",
                [
                    Resource(
                        relative_package_path = "path1.dev.js",
                        external_url = "external_path1.dev.js"
                    )
                ]
            ),
            prod = ResourcePkg(
                "dash_renderer",
                "path",
                [
                    Resource(
                        relative_package_path = "path1.prod.js",
                        external_url = "external_path1.prod.js"
                    )
                ]
            )
        ),
        dash_renderer = ResourcePkg(
            "dash_renderer",
            "path",
            [
                Resource(
                    relative_package_path = "dash-renderer/dash_renderer.min.js",
                    external_url = "https://dash_renderer.min.js"
                ),
                Resource(
                    relative_package_path = "dash-renderer/dash_renderer.dyn.js",
                    external_url = "https://dash_renderer.dyn.js",
                    dynamic = true
                ),
                Resource(
                    relative_package_path = "dash-renderer/dash_renderer.eag.js",
                    external_url = "https://dash_renderer.eag.js",
                    async = :eager
                )

            ]
        )
    )

    test_app = dash(include_assets_files = false)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 0
    @test length(app_resources.js) == 2
    @test app_resources.js[2].relative_path == "dash-renderer/dash_renderer.min.js"

    @test haskey(app_resources.files, "dash_renderer")
    @test app_resources.files["dash_renderer"].base_path == "path"
    @test "path1.prod.js" in app_resources.files["dash_renderer"].files
    @test "dash-renderer/dash_renderer.min.js" in app_resources.files["dash_renderer"].files
    @test "dash-renderer/dash_renderer.dyn.js" in app_resources.files["dash_renderer"].files
    @test "dash-renderer/dash_renderer.eag.js" in app_resources.files["dash_renderer"].files

    test_app = dash(eager_loading = true, include_assets_files = false)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 0
    @test length(app_resources.js) == 3
    @test app_resources.js[2].relative_path == "dash-renderer/dash_renderer.min.js"
    @test app_resources.js[3].relative_path == "dash-renderer/dash_renderer.eag.js"
end

@testset "application resources components" begin
    test_registry = ResourcesRegistry(
        dash_dependency = (
            dev = ResourcePkg(
                "dash_renderer",
                "path",
                [
                    Resource(
                        relative_package_path = "path1.dev.js",
                        external_url = "external_path1.dev.js"
                    )
                ]
            ),
            prod = ResourcePkg(
                "dash_renderer",
                "path",
                [
                    Resource(
                        relative_package_path = "path1.prod.js",
                        external_url = "external_path1.prod.js"
                    )
                ]
            )
        ),
        dash_renderer = ResourcePkg(
            "dash_renderer",
            "path",
            [
                Resource(
                    relative_package_path = "dash-renderer/dash_renderer.min.js",
                    external_url = "https://dash_renderer.min.js"
                ),
            ]
        )
    )

    register_package!(test_registry,
        ResourcePkg(
            "comps",
            "comps_path",
            [
                Resource(
                    relative_package_path = "comp.js"
                ),
                Resource(
                    relative_package_path = "comp.dyn.js",
                    dynamic = true
                ),
                Resource(
                    relative_package_path = "comp.css",
                    type = :css
                )
            ]
        )
    )

    test_app = dash(include_assets_files = false)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 1
    @test length(app_resources.js) == 3

    @test app_resources.js[1].relative_path == "path1.prod.js"
    @test app_resources.js[2].namespace == "comps"
    @test app_resources.js[2].relative_path == "comp.js"
    @test app_resources.js[end].relative_path == "dash-renderer/dash_renderer.min.js"

    @test haskey(app_resources.files, "comps")
    @test app_resources.files["comps"].base_path == "comps_path"
    @test "comp.js" in app_resources.files["comps"].files
    @test "comp.dyn.js" in app_resources.files["comps"].files
    @test "comp.css" in app_resources.files["comps"].files
end

@testset "application resources external" begin
    test_registry = ResourcesRegistry(
        dash_dependency = (
            dev = ResourcePkg(
                "dash_renderer",
                "path",
                [
                    Resource(
                        relative_package_path = "path1.dev.js",
                        external_url = "external_path1.dev.js"
                    )
                ]
            ),
            prod = ResourcePkg(
                "dash_renderer",
                "path",
                [
                    Resource(
                        relative_package_path = "path1.prod.js",
                        external_url = "external_path1.prod.js"
                    )
                ]
            )
        ),
        dash_renderer = ResourcePkg(
            "dash_renderer",
            "path",
            [
                Resource(
                    relative_package_path = "dash-renderer/dash_renderer.min.js",
                    external_url = "https://dash_renderer.min.js"
                ),
            ]
        )
    )

    external_css = ["test.css", Dict("rel"=>"stylesheet", "href"=>"test.css")]
    external_js = ["test.js", Dict("crossorigin"=>"anonymous", "src"=>"test.js")]
    test_app = dash(external_stylesheets = external_css, include_assets_files = false, external_scripts = external_js)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 2
    @test length(app_resources.js) == 4

    @test app_resources.js[1].relative_path == "path1.prod.js"
    @test app_resources.js[end].relative_path == "dash-renderer/dash_renderer.min.js"

    @test app_resources.css[1] isa AppExternalResource
    @test app_resources.css[1].url == "test.css"

    @test app_resources.css[2] isa AppCustomResource
    @test app_resources.css[2].params["href"] == "test.css"
    @test app_resources.css[2].params["rel"] == "stylesheet"

    @test app_resources.js[2] isa AppExternalResource
    @test app_resources.js[2].url == "test.js"

    @test app_resources.js[3] isa AppCustomResource
    @test app_resources.js[3].params["src"] == "test.js"
    @test app_resources.js[3].params["crossorigin"] == "anonymous"

end
@testset "application resources assets" begin
    test_registry = ResourcesRegistry(
        dash_dependency = (
            dev = ResourcePkg(
                "dash_renderer",
                "path",
                [
                ]
            ),
            prod = ResourcePkg(
                "dash_renderer",
                "path",
                [
                ]
            )
        ),
        dash_renderer = ResourcePkg(
            "dash_renderer",
            "path",
            [
            ]
        )
    )

    test_app = dash(assets_folder = "assets_include")

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 1
    @test length(app_resources.js) == 4
    @test app_resources.favicon == "favicon.ico"
    @test all(app_resources.css) do r
        r isa AppAssetResource
    end
    @test app_resources.css[1].path == "test.css"

    @test all(app_resources.js) do r
        r isa AppAssetResource
    end

    @test any(app_resources.js) do r
        r.path == "test1.js"
    end
    @test any(app_resources.js) do r
        r.path == "test2.js"
    end
    @test any(app_resources.js) do r
        r.path == "sub/sub1.js"
    end
    @test any(app_resources.js) do r
        r.path == "sub/sub2.js"
    end

    test_app = dash(assets_folder = "assets_include", include_assets_files = false)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 0
    @test length(app_resources.js) == 0
    @test app_resources.favicon ===  nothing

    test_app = dash(assets_folder = "assets_include", serve_locally = false)

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 1
    @test length(app_resources.js) == 4
    @test app_resources.favicon == "favicon.ico"
    @test all(app_resources.css) do r
        r isa AppAssetResource
    end
    @test app_resources.css[1].path == "test.css"

    @test all(app_resources.js) do r
        r isa AppAssetResource
    end

    @test any(app_resources.js) do r
        r.path == "test1.js"
    end
    @test any(app_resources.js) do r
        r.path == "test2.js"
    end
    @test any(app_resources.js) do r
        r.path == "sub/sub1.js"
    end
    @test any(app_resources.js) do r
        r.path == "sub/sub2.js"
    end

    test_app = dash(assets_folder = "assets_include", serve_locally = false, assets_external_path = "http://ext/")

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 1
    @test length(app_resources.js) == 4
    @test app_resources.favicon == "favicon.ico"

    @test all(app_resources.css) do r
        r isa AppExternalResource
    end

    @test app_resources.css[1].url == "http://ext/test.css"

    @test all(app_resources.js) do r
        r isa AppExternalResource
    end

    @test any(app_resources.js) do r
        r.url == "http://ext/test1.js"
    end
    @test any(app_resources.js) do r
        r.url == "http://ext/test2.js"
    end
    @test any(app_resources.js) do r
        r.url == "http://ext/sub/sub1.js"
    end
    @test any(app_resources.js) do r
        r.url == "http://ext/sub/sub2.js"
    end

    test_app = dash(assets_folder = "assets_include", assets_ignore = ".*1")

    app_resources = ApplicationResources(test_app, test_registry)
    @test length(app_resources.css) == 1
    @test length(app_resources.js) == 2
    @test app_resources.favicon == "favicon.ico"
    @test all(app_resources.css) do r
        r isa AppAssetResource
    end
    @test app_resources.css[1].path == "test.css"

    @test all(app_resources.js) do r
        r isa AppAssetResource
    end

    @test !any(app_resources.js) do r
        r.path == "test1.js"
    end
    @test any(app_resources.js) do r
        r.path == "test2.js"
    end
    @test !any(app_resources.js) do r
        r.path == "sub/sub1.js"
    end
    @test any(app_resources.js) do r
        r.path == "sub/sub2.js"
    end
end


text:
@testset "broadcast" begin
    i = Indices([1,2,3,4,5])

    @test isequal((i .+ 1)::Dictionary, dictionary([1=>2, 2=>3, 3=>4, 4=>5, 5=>6]))
    @test isequal(Base.Broadcast.broadcasted(+, i, 1)::BroadcastedDictionary, dictionary([1=>2, 2=>3, 3=>4, 4=>5, 5=>6]))

    @test isequal((i .+ i)::Dictionary, dictionary([1=>2, 2=>4, 3=>6, 4=>8, 5=>10]))
    @test isequal(Base.Broadcast.broadcasted(+, i, i)::BroadcastedDictionary, dictionary([1=>2, 2=>4, 3=>6, 4=>8, 5=>10]))

    @test isequal((i .+ copy(i))::Dictionary, dictionary([1=>2, 2=>4, 3=>6, 4=>8, 5=>10]))
    @test isequal(Base.Broadcast.broadcasted(+, i, copy(i))::BroadcastedDictionary, dictionary([1=>2, 2=>4, 3=>6, 4=>8, 5=>10]))

    @test_throws IndexError Indices([1,2]) .+ Indices([2,3])

    d = i .+ 1

    @test isequal((d .+ 1)::Dictionary, dictionary([1=>3, 2=>4, 3=>5, 4=>6, 5=>7]))
    @test isequal(Base.Broadcast.broadcasted(+, d, 1)::BroadcastedDictionary, dictionary([1=>3, 2=>4, 3=>5, 4=>6, 5=>7]))

    @test isequal((d .+ d)::Dictionary, dictionary([1=>4, 2=>6, 3=>8, 4=>10, 5=>12]))
    @test isequal(Base.Broadcast.broadcasted(+, d, d)::BroadcastedDictionary, dictionary([1=>4, 2=>6, 3=>8, 4=>10, 5=>12]))

    @test isequal((d .+ copy(d))::Dictionary, dictionary([1=>4, 2=>6, 3=>8, 4=>10, 5=>12]))
    @test isequal(Base.Broadcast.broadcasted(+, d, copy(d))::BroadcastedDictionary, dictionary([1=>4, 2=>6, 3=>8, 4=>10, 5=>12]))

    d2 = similar(d)
    d2 .= d .+ d
    @test isequal(d2, dictionary([1=>4, 2=>6, 3=>8, 4=>10, 5=>12]))

    @test_throws IndexError Dictionary([1,2],[1,2]) .+ Dictionary([2,3],[2,3])
end


text:
"""
    free_energy!(model, config)

Compute the free energy of the SCFT polymer model. The chemical potential is also returned.

Only applicable to two-chain two-component system, e.g. AB/A, AB3/A etc.
"""
function free_energy!(model, config)
    cellopt = uppercase(config["Algorithm_Cell_Optimization"]["algorithm"])
    if cellopt == "SINGLE"
        F = run!(model, config)
    elseif cellopt == "BRENT" || cellopt == "AUTO"
        res = optimize_cell!(model, config)
        config["Algorithm_Cell_Optimization"]["cell_guess"][1] = Optim.minimizer(res)
        # F = minimum(res)
    end

    return Polyorder.H(model)+Polyorder.Hi(model), Polyorder.mu(model)
end

"""
    grand_potential(ϕ₁, ϕ₂, ϕ₀, model1, config1; model2=nothing, config2=nothing)

Compute the grand potential and its first order derivatives. Note the same polymer model is assumed is for model1 and model2 but for different phases and chain volume fractions. If model2 and config2 are also specified. Theire consistency with model1 and config1 is not checked! Please make sure they describe the same polymer model!

Only applicable to two-chain two-component system, e.g. AB/A, AB3/A etc.

This method is deprecated and replaced by PhaseDiagram.jl: gibbs_free_energy.
"""
function grand_potential(ϕ₁, ϕ₂, ϕ₀, model1, config1; model2=nothing, config2=nothing)
    if model2 === nothing
        model2 = model1
    end
    if config2 === nothing
        config2 = config1
    end
    v₁ = (ϕ₀ - ϕ₂) / (ϕ₁ - ϕ₂)
    v₂ = 1 - v₁
    config1["Model"]["phi"][1] = 1 - ϕ₁
    config1["Model"]["phi"][2] = ϕ₁
    F₁, μ₁ = free_energy!(model1, config1)
    config2["Model"]["phi"][1] = 1 - ϕ₂
    config2["Model"]["phi"][2] = ϕ₂
    F₂, μ₂ = free_energy!(model2, config2)
    G = v₁ * F₁ + v₂ * F₂
    fs = config1["Model"]["f"]
    α = fs[1] * fs[end]
    dG1 = v₁ * (F₂ - F₁) / (ϕ₁ - ϕ₂) + v₁ * μ₁ / α
    dG2 = v₂ * (F₂ - F₁) / (ϕ₁ - ϕ₂) + v₂ * μ₂ / α
    return G, dG1, dG2
end


text:
function isfair(G::GraphGame, L::NTuple{2,Int})
    isfair(G, Edge(L))
end
function isfair(G::GraphGame, L::Edge)

    g′ = copy(G.g)
    rem_edge!(g′, L)

    G′ = GraphGame(G.N, G.v, g′)

    i, j = L.src, L.dst

    f′ = fₛ(G′)
    f = fₛ(G)

    return f[i] - f′[i] ≈ f[j] - f′[j]

end


text:
module NPC

using ..Ahorn, Maple

const placements = Ahorn.PlacementDict(
    "NPC" => Ahorn.EntityPlacement(
        Maple.NPC
    )
)

Ahorn.editingOptions(entity::Maple.NPC) = Dict{String, Any}(
    "npc" => Maple.npc_npcs
)

npcSprites = Dict{String, String}(
    "granny" => "characters/oldlady/idle00",
    "theo" => "characters/theo/theo00",
    "oshiro" => "characters/oshiro/oshiro24",
    "evil" => "characters/badeline/sleep00",
    "badeline" => "characters/badeline/sleep00",
)

fallbackSprite = "characters/oldlady/idle00"

function getTexture(entity::Maple.NPC)
    npcName = get(entity.data, "npc", "granny_00_house")
    name = lowercase(split(npcName, "_")[1])

    return get(npcSprites, name, false)
end

function Ahorn.selection(entity::Maple.NPC)
    x, y = Ahorn.position(entity)
    sprite = getTexture(entity)

    if isa(sprite, String)
        return Ahorn.getSpriteRectangle(sprite, x, y, jx=0.5, jy=1.0)
    end

    return Ahorn.getSpriteRectangle(fallbackSprite, x, y, jx=0.5, jy=1.0)
end

function Ahorn.render(ctx::Ahorn.Cairo.CairoContext, entity::Maple.NPC, room::Maple.Room)
    texture = getTexture(entity)

    if isa(texture, String)
        Ahorn.drawSprite(ctx, texture, 0, 0, jx=0.5, jy=1.0)
    end

    return true
end

end


text:
__precompile__(true)

module Giza

export
    GIZA_LS_SOLID,
    GIZA_LS_SHORT_DASH,
    GIZA_LS_DASHED,
    GIZA_LS_LONG_DASH,
    GIZA_LS_DOT,
    GIZA_LS_DOTTED,
    GIZA_LS_DASH_DOT,
    GIZA_LS_DASH_DOT_DOT_DOT,

    GIZA_NUMBER_FORMAT_AUTO,
    GIZA_NUMBER_FORMAT_DEC,
    GIZA_NUMBER_FORMAT_EXP,

    GIZA_UNITS_NORMALISED,
    GIZA_UNITS_INCHES,
    GIZA_UNITS_MM,
    GIZA_UNITS_PIXELS,
    GIZA_UNITS_WORLD,
    GIZA_UNITS_DEVICE,

    GIZA_BACKGROUND,
    GIZA_BACKGROUND_COLOUR,
    GIZA_BACKGROUND_COLOR,
    GIZA_FOREGROUND,
    GIZA_FOREGROUND_COLOUR,
    GIZA_FOREGROUND_COLOR,

    GIZA_COLOUR_INDEX_MIN,
    GIZA_COLOUR_INDEX_MAX,

    GIZA_COLOUR_PALETTE_DEFAULT,
    GIZA_COLOR_PALETTE_DEFAULT,
    GIZA_COLOUR_PALETTE_PGPLOT,
    GIZA_COLOR_PALETTE_PGPLOT,

    GIZA_FILL_SOLID,
    GIZA_FILL_HOLLOW,
    GIZA_FILL_HATCH,
    GIZA_FILL_CROSSHATCH,

    GIZA_MAX_FILL_STYLES,

    GIZA_EXTEND_NONE,
    GIZA_EXTEND_REPEAT,
    GIZA_EXTEND_REFLECT,
    GIZA_EXTEND_PAD,

    GIZA_BAND_NONE,
    GIZA_BAND_LINE,
    GIZA_BAND_RECTANGLE,
    GIZA_BAND_HORZLINES,
    GIZA_BAND_VERTLINES,
    GIZA_BAND_HORZLINE,
    GIZA_BAND_VERTLINE,
    GIZA_BAND_CROSSHAIR,
    GIZA_BAND_CIRCLE,

    GIZA_LEFT_CLICK,
    GIZA_MIDDLE_CLICK,
    GIZA_RIGHT_CLICK,
    GIZA_SHIFT_CLICK,
    GIZA_SCROLL_UP,
    GIZA_SCROLL_DOWN,
    GIZA_SCROLL_LEFT,
    GIZA_SCROLL_RIGHT,
    GIZA_OTHER_CLICK

if isfile(joinpath(@__DIR__, "..", "deps", "deps.jl"))
    include(joinpath("..", "deps", "deps.jl"))
else
    if false
        # For official package.
        error("Giza not properly installed.  Please run Pkg.build(\"Giza\")")
    else
        # For un-official package.
        error(string("Giza not properly installed.  Please run:\n\n",
                     "    cd \"", normpath(@__DIR__, "..", "deps"), "\"\n",
                     "    make GIZA_JL_DLL=\"\$PATH_TO_GIZA_DLL\"\n"))
    end
end

_array(::Type{T}, A::DenseArray{T,N}) where {T<:Real,N} = A
_array(::Type{T}, A::AbstractArray{<:Real,N}) where {T<:Real,N} =
    convert(Array{T,N}, A)

_float(::Type{Cfloat}) = Cfloat
_float(::Type{<:Real}) = Cdouble
_float(::Type{Cfloat}, ::Type{Cfloat}) = Cfloat
_float(::Type{<:Real}, ::Type{<:Real}) = Cdouble
_float(::Type{Cfloat}, ::Type{Cfloat}, ::Type{Cfloat}) = Cfloat
_float(::Type{<:Real}, ::Type{<:Real}, ::Type{<:Real}) = Cdouble

const GIZA_LS_SOLID            = Cint(1)
const GIZA_LS_SHORT_DASH       = Cint(2)
const GIZA_LS_DASHED           = Cint(2)
const GIZA_LS_LONG_DASH        = Cint(3)
const GIZA_LS_DOT              = Cint(4)
const GIZA_LS_DOTTED           = Cint(4)
const GIZA_LS_DASH_DOT         = Cint(5)
const GIZA_LS_DASH_DOT_DOT_DOT = Cint(6)

const GIZA_NUMBER_FORMAT_AUTO  = Cint(0)
const GIZA_NUMBER_FORMAT_DEC   = Cint(1)
const GIZA_NUMBER_FORMAT_EXP   = Cint(2)

const GIZA_UNITS_NORMALISED = Cint(0)
const GIZA_UNITS_INCHES     = Cint(1)
const GIZA_UNITS_MM         = Cint(2)
const GIZA_UNITS_PIXELS     = Cint(3)
const GIZA_UNITS_WORLD      = Cint(4)
const GIZA_UNITS_DEVICE     = Cint(5)

const GIZA_BACKGROUND        = Cint(0)
const GIZA_BACKGROUND_COLOUR = GIZA_BACKGROUND
const GIZA_BACKGROUND_COLOR  = GIZA_BACKGROUND
const GIZA_FOREGROUND        = Cint(1)
const GIZA_FOREGROUND_COLOUR = GIZA_FOREGROUND
const GIZA_FOREGROUND_COLOR  = GIZA_FOREGROUND

const GIZA_COLOUR_INDEX_MIN  =   0
const GIZA_COLOUR_INDEX_MAX  = 271

const GIZA_COLOUR_PALETTE_DEFAULT = Cint(0)
const GIZA_COLOR_PALETTE_DEFAULT  = GIZA_COLOUR_PALETTE_DEFAULT
const GIZA_COLOUR_PALETTE_PGPLOT  = Cint(1)
const GIZA_COLOR_PALETTE_PGPLOT   = GIZA_COLOUR_PALETTE_PGPLOT

const GIZA_FILL_SOLID      = Cint(1)
const GIZA_FILL_HOLLOW     = Cint(2)
const GIZA_FILL_HATCH      = Cint(3)
const GIZA_FILL_CROSSHATCH = Cint(4)

const GIZA_MAX_FILL_STYLES = 4

const GIZA_EXTEND_NONE    = Cint(0)
const GIZA_EXTEND_REPEAT  = Cint(1)
const GIZA_EXTEND_REFLECT = Cint(2)
const GIZA_EXTEND_PAD     = Cint(3)

function arrow(x1::Real, y1::Real, x2::Real, y2::Real)
    ccall((:giza_arrow, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble),
          x1, y1, x2, y2)
end

function set_arrow_style(fillstyle::Integer, angle::Real, cutback::Real)
     ccall((:giza_set_arrow_style, _GIZALIB), Cvoid,
           (Cint, Cdouble, Cdouble),
           fillstyle, angle, cutback)
end

function get_arrow_style()
    fillstyle = Ref{Cint}()
    angle = Ref{Cdouble}()
    fillstyle = Ref{Cdouble}()
     ccall((:giza_get_arrow_style, _GIZALIB), Cvoid,
           (Ptr{Cint}, Ptr{Cdouble}, Ptr{Cdouble}),
           fillstyle, angle, cutback)
    return (fillstyle[], angle[], cutback[])
end

function annotate(side::AbstractString,
                  offset::Real,
                  coord::Real,
                  just::Real,
                  text::AbstractString)
    ccall((:giza_annotate, _GIZALIB), Cvoid,
          (Cstring, Cdouble, Cdouble, Cdouble, Cstring),
          side, offset, coord, just, text)
end

function label(xlab::AbstractString,
               ylab::AbstractString,
               title::AbstractString)
    ccall((:giza_label, _GIZALIB), Cvoid,
          (Cstring, Cstring, Cstring),
          xlab, ylab, title)
end

const GIZA_BAND_NONE      = Cint(0)
const GIZA_BAND_LINE      = Cint(1)
const GIZA_BAND_RECTANGLE = Cint(2)
const GIZA_BAND_HORZLINES = Cint(3)
const GIZA_BAND_VERTLINES = Cint(4)
const GIZA_BAND_HORZLINE  = Cint(5)
const GIZA_BAND_VERTLINE  = Cint(6)
const GIZA_BAND_CROSSHAIR = Cint(7)
const GIZA_BAND_CIRCLE    = Cint(8)

const GIZA_LEFT_CLICK     = 'A'
const GIZA_MIDDLE_CLICK   = 'D'
const GIZA_RIGHT_CLICK    = 'X'
const GIZA_SHIFT_CLICK    = Char(15)
const GIZA_SCROLL_UP      = Char(21)
const GIZA_SCROLL_DOWN    = Char( 4)
const GIZA_SCROLL_LEFT    = Char(12)
const GIZA_SCROLL_RIGHT   = Char(18)
const GIZA_OTHER_CLICK    = 'X'

function band(mode::Integer, move::Bool, xanc::Real, yanc::Real)
    x = Ref{Cdouble}()
    y = Ref{Cdouble}()
    c = Ref{Cchar}()
    code = ccall((:giza_band, _GIZALIB), Cint,
                 (Cint, Cint, Cdouble, Cdouble,
                  Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cchar}),
                 mode, move, xanc, yanc, x, y, c)
    if code != 0
        error("giza_band failed with error $code")
    end
    return (x[], y[], Char(c[]))
end

function get_band_style()
    ls = Ref{Cint}()
    lw = Ref{Cdouble}()
    ccall((:giza_get_band_style, _GIZALIB), Cvoid,
          (Ptr{Cint}, Ptr{Cdouble}), ls, lw)
    return (ls[], lw[])
end

function set_band_style(ls::Integer, lw::Real)
    ccall((:giza_set_band_style, _GIZALIB), Cvoid,
          (Cint, Cdouble), ls, lw)
end

for f in (:box, :box_time)
    @eval function $f(xopt::AbstractString, xtick::Real, nxsub::Integer,
	              yopt::AbstractString, ytick::Real, nysub::Integer)
        ccall(($(string("giza_",f)), _GIZALIB), Cvoid,
              (Cstring, Cdouble, Cint, Cstring, Cdouble, Cint),
              xopt, xtick, nxsub, yopt, ytick, nysub)
    end
end

function get_buffering()
    buf = Ref{Cint}()
    ccall((:giza_get_buffering, _GIZALIB), Cvoid, (Ptr{Cint},), buf)
    return (buf[] != 0)
end

function circle(x::Real, y::Real, r::Real)
    ccall((:giza_circle, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble),
          x, y, r)
end

set_character_height(h::Real) =
    ccall((:giza_set_character_height, _GIZALIB), Cvoid, (Cdouble,), h)

function get_character_height()
    h = Ref{Cdouble}()
    ccall((:giza_get_character_height, _GIZALIB), Cvoid,
          (Ptr{Cdouble},), h)
    return h[]
end

function get_character_size(units::Integer)
    x = Ref{Cdouble}()
    y = Ref{Cdouble}()
    ccall((:giza_get_character_size, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}), units, x, y)
    return (x[], y[])
end

function set_clipping(clip::Bool)
    ccall((:giza_set_clipping, _GIZALIB), Cvoid, (Cint,), clip)
end

function get_clipping()
    clip = Ref{Cint}()
    ccall((:giza_get_clipping, _GIZALIB), Cvoid, (Ptr{Cint},), clip)
    return (clip[] != 0)
end

color_bar(args...) = set_colour_bar(args...)

function colour_bar(side::AbstractString, disp::Real, width::Real,
                    vmin::Real, vmax::Real, label::AbstractString)
    ccall((:giza_colour_bar, _GIZALIB), Cvoid,
          (Cstring, Cdouble, Cdouble, Cdouble, Cdouble, Cstring),
          side, disp, width, vmin, vmax, label)
end

function set_colour_index(ci::Integer)
    ccall((:giza_set_colour_index, _GIZALIB), Cvoid, (Cint,), ci)
end

function get_colour_index()
    ci = Ref{Cint}()
    ccall((:giza_get_colour_index, _GIZALIB), Cvoid, (Ptr{Cint},), ci)
    return ci[]
end

function set_colour_index_range(cimin::Integer, cimax::Integer)
    ccall((:giza_set_colour_index_range, _GIZALIB), Cvoid,
          (Cint, Cint), cimin, cimax)
end

function get_colour_index_range()
    cimin = Ref{Cint}()
    cimax = Ref{Cint}()
    ccall((:giza_get_colour_index_range, _GIZALIB), Cvoid,
          (Ptr{Cint}, Ptr{Cint}), cimin, cimax)
    return (cimin[], cimax[])
end

function set_colour_palette(pal::Integer)
    ccall((:giza_set_colour_palette, _GIZALIB), Cvoid, (Cint,), pal)
end

function set_colour_representation(ci::Integer,
                                   red::Real,
                                   green::Real,
                                   blue::Real)
    ccall((:giza_set_colour_representation, _GIZALIB), Cvoid,
          (Cint, Cdouble, Cdouble, Cdouble), ci, red, green, blue)
end

function get_colour_representation(ci::Integer)
    red   = Ref{Cdouble}()
    green = Ref{Cdouble}()
    blue  = Ref{Cdouble}()
    ccall((:giza_get_colour_representation, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}),
          ci, red, green, blue)
    return (red[], green[], blue[])
end

function set_colour_representation_alpha(ci::Integer,
                                         red::Real,
                                         green::Real,
                                         blue::Real,
                                         alpha::Real)
    ccall((:giza_set_colour_representation_alpha, _GIZALIB), Cvoid,
          (Cint, Cdouble, Cdouble, Cdouble, Cdouble),
          ci, red, green, blue, alpha)
end

function get_colour_representation_alpha(ci::Integer)
    red   = Ref{Cdouble}()
    green = Ref{Cdouble}()
    blue  = Ref{Cdouble}()
    alpha = Ref{Cdouble}()
    ccall((:giza_get_colour_representation_alpha, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}),
          ci, red, green, blue, alpha)
    return (red[], green[], blue[], alpha[])
end

function set_colour_representation_rgb(ci::Integer,
                                       red::Integer,
                                       green::Integer,
                                       blue::Integer)
    ccall((:giza_set_colour_representation_rgb, _GIZALIB), Cvoid,
          (Cint, Cint, Cint, Cint), ci, red, green, blue)
end

function get_colour_representation_rgb(ci::Integer)
    red   = Ref{Cint}()
    green = Ref{Cint}()
    blue  = Ref{Cint}()
   ccall((:giza_get_colour_representation_rgb, _GIZALIB), Cvoid,
          (Cint, Ptr{Cint}, Ptr{Cint}, Ptr{Cint}),
          ci, red, green, blue)
    return (red[], green[], blue[])
end

function set_colour_representation_rgba(ci::Integer,
                                        red::Integer,
                                        green::Integer,
                                        blue::Integer,
                                        alpha::Real)
    ccall((:giza_set_colour_representation_rgba, _GIZALIB), Cvoid,
          (Cint, Cint, Cint, Cint, Cdouble), ci, red, green, blue, alpha)
end

function get_colour_representation_rgba(ci::Integer)
    red   = Ref{Cint}()
    green = Ref{Cint}()
    blue  = Ref{Cint}()
    alpha = Ref{Cdouble}()
    ccall((:giza_get_colour_representation_rgba, _GIZALIB), Cvoid,
          (Cint, Ptr{Cint}, Ptr{Cint}, Ptr{Cint}, Ptr{Cdouble}),
          ci, red, green, blue, alpha)
    return (red[], green[], blue[], alpha[])
end

function set_colour_representation_hls(ci::Integer,
                                       hue::Real,
                                       lightness::Real,
                                       saturation::Real)
    ccall((:giza_set_colour_representation_hls, _GIZALIB), Cvoid,
          (Cint, Cdouble, Cdouble, Cdouble),
          ci, hue, lightness, saturation)
end

function set_colour_table(pos::AbstractVector{<:Real},
                          red::AbstractVector{<:Real},
                          green::AbstractVector{<:Real},
                          blue::AbstractVector{<:Real},
                          contrast::Real,
                          brightness::Real)
    n = length(pos)
    @assert length(red) == n
    @assert length(green) == n
    @assert length(blue) == n
    code = ccall((:giza_set_colour_table, _GIZALIB), Cint,
                 (Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble},
                  Cint, Cdouble, Cdouble),
                 _array(Cdouble, pos), _array(Cdouble, red),
                 _array(Cdouble, green), _array(Cdouble, blue),
                 n, contrast, brightness)
    if code != 0
        error("giza_set_colour_table failed with error $code")
    end
end

function set_colour_table_gray()
    code = ccall((:giza_set_colour_table_gray, _GIZALIB), Cint, ())
    if code != 0
        error("giza_set_colour_table_gray failed with error $code")
    end
end

function rgb_from_table(pos::Real)
    red   = Ref{Cdouble}()
    green = Ref{Cdouble}()
    blue  = Ref{Cdouble}()
    ccall((:giza_rgb_from_table, _GIZALIB), Cvoid,
          (Cdouble, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}),
          pos, red, green, blue)
    return (red[], green[], blue[])
end

function contour(arr::AbstractMatrix{<:Real},
                 cont::AbstractVector{<:Real},
                 aff::AbstractArray{<:Real,N}) where {N}
    @assert length(aff) == 6
    xdim, ydim = size(arr)
    ccall((:giza_contour, _GIZALIB), Cvoid,
          (Cint, Cint, Ptr{Cdouble}, Cint, Cint, Cint, Cint, Cint,
	   Ptr{Cdouble}, Ptr{Cdouble}),
          xdim, ydim, _array(Cdouble, arr), 0, xdim - 1, 0, ydim - 1,
          length(cont), _array(Cdouble, cont), _array(Cdouble, affine))
end

function get_current_point()
    x = Ref{Cdouble}()
    y = Ref{Cdouble}()
    ccall((:giza_get_current_point, _GIZALIB), Cvoid,
          (Ptr{Cdouble}, Ptr{Cdouble}), x, y)
    return (x[], y[])
end

function open_device(devname::AbstractString, prefix::AbstractString)
    id = ccall((:giza_open_device, _GIZALIB), Cint,
               (Cstring, Cstring), devname, prefix)
    if id <= 0
        error("giza_open_device failed")
    end
    return id
end

function open_device_size(devname::AbstractString, prefix::AbstractString,
                          width::Real, height::Real, units::Integer)
    code = ccall((:giza_open_device_size, _GIZALIB), Cint,
          (Cstring, Cstring, Cdouble, Cdouble, Cint),
          devname, prefix, width, height, units)
    if code != 0
        error("giza_open_device falied with error $code")
    end
end

function select_device(dev::Integer)
    ccall((:giza_select_device, _GIZALIB), Cvoid, (Cint,), dev)
end

function get_device_id()
    dev = Ref{Cint}()
    ccall((:giza_get_device_id, _GIZALIB), Cvoid, (Ptr{Cint},), dev)
    return dev[]
end

function query_device(key::AbstractString)
    buf = Array{UInt8}(undef, 3000) # FIXME: the function below is badly coded
    ccall((:giza_query_device, _GIZALIB), Cint,
          (Cstring, Ptr{UInt8}), key, buf)
    buf[end] = 0
    return unsafe_string(pointer(buf))
end

function device_has_cursor()
    return (ccall((:giza_device_has_cursor, _GIZALIB), Cint, ()) != 0)
end

function get_key_press()
    x = Ref{Cdouble}()
    y = Ref{Cdouble}()
    c = Ref{Cchar}()
    code = ccall((:giza_get_key_press, _GIZALIB), Cint,
                 (Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cchar}),
                 x, y, c)
    if code != 0
        error("giza_get_key_press failed with error $code")
    end
    return (x[], y[], Char(c[]))
end

function draw(x::Real, y::Real)
    ccall((:giza_draw, _GIZALIB), Cvoid,
          (Cdouble, Cdouble), x, y)
end

function set_environment(xmin::Real, xmax::Real,
                         ymin::Real, ymax::Real,
                         just::Integer, axis::Integer)
    ccall((:giza_set_environment, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble, Cint, Cint),
          xmin, xmax, ymin, ymax, just, axis)
end

function error_bars(dir::Integer,
                    xpts::AbstractVector{Tx},
                    ypts::AbstractVector{Ty},
                    errs::AbstractVector{Te},
                    term::Real) where {Tx<:Real, Ty<:Real, Te<:Real}
    @assert length(xpts) == length(ypts) == length(errs)
    n = length(xpts)
    T = _float(Tx, Ty, Te)
    _error_bars(dir, n,
                _array(T, xpts),
                _array(T, ypts),
                _array(T, errs), term)
end

function _error_bars(dir::Integer,
                     n::Integer,
                     xpts::DenseVector{Cdouble},
                     ypts::DenseVector{Cdouble},
                     errs::DenseVector{Cdouble},
                     term::Real)
    ccall((:giza_error_bars, _GIZALIB), Cvoid,
          (Cint, Cint, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}, Cdouble),
          dir, n, xpts, ypts, errs, term)
end

function _error_bars(dir::Integer,
                     n::Integer,
                     xpts::DenseVector{Cfloat},
                     ypts::DenseVector{Cfloat},
                     errs::DenseVector{Cfloat},
                     term::Cfloat)
    ccall((:giza_error_bars_float, _GIZALIB), Cvoid,
          (Cint, Cint, Ptr{Cfloat}, Ptr{Cfloat}, Ptr{Cfloat}, Cfloat),
          dir, n, xpts, ypts, errs, term)
end

function error_bars_vert(dir::Integer,
                         xpts::AbstractVector{Tx},
                         ypts1::AbstractVector{Ty1},
                         ypts2::AbstractVector{Ty2},
                         term::Real) where {Tx<:Real, Ty1<:Real, Ty2<:Real}
    @assert length(xpts) == length(ypts1) == length(ypts2)
    n = length(xpts)
    T = _float(Tx, Ty1, Ty2)
    _error_bars_vert(dir, n,
                     _array(T, xpts),
                     _array(T, ypts1),
                     _array(T, ypts2), term)
end

function _error_bars_vert(n::Integer,
                          xpts::DenseVector{Cdouble},
                          ypts1::DenseVector{Cdouble},
                          ypts2::DenseVector{Cdouble},
                          term::Real)
    ccall((:giza_error_bars_vert, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}, Cdouble),
          n, xpts, ypts1, ypts2, term)
end

function _error_bars_vert(n::Integer,
                          xpts::DenseVector{Cfloat},
                          ypts1::DenseVector{Cfloat},
                          ypts2::DenseVector{Cfloat},
                          term::Cfloat)
    ccall((:giza_error_bars_vert_float, _GIZALIB), Cvoid,
          (Cint, Ptr{Cfloat}, Ptr{Cfloat}, Ptr{Cfloat}, Cfloat),
          n, xpts, ypts1, ypts2, term)
end

function error_bars_hori(dir::Integer,
                         xpts1::AbstractVector{Tx1},
                         xpts2::AbstractVector{Tx2},
                         ypts::AbstractVector{Ty},
                         term::Real) where {Tx1<:Real, Tx2<:Real, Ty<:Real}
    @assert length(xpts1) == length(xpts2) == length(ypts)
    n = length(xpts)
    T = _float(Tx1, Tx2, Ty)
    _error_bars_hori(dir, n,
                     _array(T, xpts1),
                     _array(T, xpts2),
                     _array(T, ypts), term)
end

function _error_bars_hori(n::Integer,
                          xpts1::DenseVector{Cdouble},
                          xpts2::DenseVector{Cdouble},
                          ypts::DenseVector{Cdouble},
                          term::Real)
    ccall((:giza_error_bars_hori, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}, Cdouble),
          n, xpts1, xpts2, ypts, term)
end

function _error_bars_hori(n::Integer,
                          xpts1::DenseVector{Cfloat},
                          xpts2::DenseVector{Cfloat},
                          ypts::DenseVector{Cfloat},
                          term::Cfloat)
    ccall((:giza_error_bars_hori_float, _GIZALIB), Cvoid,
          (Cint, Ptr{Cfloat}, Ptr{Cfloat}, Ptr{Cfloat}, Cfloat),
          n, xpts1, xpts2, ypts, term)
end

function set_fill(fs::Integer)
    ccall((:giza_set_fill, _GIZALIB), Cvoid, (Cint,), fs)
end

function get_fill()
    fs = Ref{Cint}()
    ccall((:giza_get_fill, _GIZALIB), Cvoid, (Ptr{Cint},), fs)
    return fs[]
end

function format_number(mantissa::Integer, power::Integer, form::Integer)
    buf = Array{UInt8}(undef, 102)
    ccall((:giza_format_number, _GIZALIB), Cvoid,
          (Cint, Cint, Cint, Ptr{UInt8}),
          mantissa, power, form, buf)
    buf[end] = 0
    return unsafe_string(pointer(buf))
end

# FIXME: Implement bindings for these?
# void giza_function_x(double (*func) (double *x), int n, double xmin,
#                     double xmax, int flag);
# void giza_function_x_float(float (*func) (float *x), int n,
#                            float xmin, float xmax, int flag);
# void giza_function_y(double (*func) (double *y), int n, double ymin,
#                      double ymax, int flag);
# void giza_function_y_float(float (*func) (float *y), int n,
#                            float ymin, float ymax, int flag);
# void giza_function_t(double (*funcx) (double *t),
#                      double (*funcy) (double *t), int n, double ymin,
#                      double ymax, int flag);
# void giza_function_t_float(float (*funcx) (float *t),
#                            float (*funcy) (float *t), int n,
#                            float ymin, float ymax, int flag);

function get_surface_size()
    x1 = Ref{Cdouble}()
    x2 = Ref{Cdouble}()
    y1 = Ref{Cdouble}()
    y2 = Ref{Cdouble}()
    ccall((:giza_get_surface_size, _GIZALIB), Cvoid,
          (Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}),
          x1, x2, y1, y2)
    return (x1[], x2[], y1[], y2[])
end

function set_hatching_style(angle::Real, spacing::Real, phase::Real)
    ccall((:giza_set_hatching_style, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble),
          angle, spacing, phase)
end

function get_hatching_style()
    angle   = Ref{Cdouble}()
    spacing = Ref{Cdouble}()
    phase   = Ref{Cdouble}()
    ccall((:giza_get_hatching_style, _GIZALIB), Cvoid,
          (Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}),
          angle, spacing, phase)
    return (angle[], spacing[], phase[])
end

function histogram(dat::AbstractArray{<:Real},
                   vmin::Real, vmax::Real,
                   nbin::Integer, flag::Integer)
    n = length(dat)
    ccall((:giza_histogram, _GIZALIB), Cvoid,
           (Cint, Ptr{Cdouble}, Cdouble, Cdouble, Cint, Cint),
          n, _array(Cdouble, dat), vmin, vmax, nbin, flag)
end

function histogram_binned(x::AbstractVector{<:Real},
                          y::AbstractVector{<:Real},
                          center::Bool)
    @assert length(x) == length(y)
    n = length(x)
    ccall((:giza_histogram_binned, _GIZALIB), Cvoid,
           (Cint, Ptr{Cdouble}, Ptr{Cdouble}, Cint),
          n, _array(Cdouble, x), _array(Cdouble, y), center)
end

function line(xpts::AbstractVector{Tx},
              ypts::AbstractVector{Ty}) where {Tx<:Real, Ty<:Real}
    @assert length(xpts) == length(xpts)
    T = _float(Tx, Ty)
    n = length(xpts)
    _line(n, _array(T, xpts), _array(T, ypts))
end

_line(n::Integer, xpts::DenseVector{Cdouble}, ypts::DenseVector{Cdouble}) =
    ccall((:giza_line, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}), n, xpts, ypts)

_line(n::Integer, xpts::DenseVector{Cfloat}, ypts::DenseVector{Cfloat}) =
    ccall((:giza_line_float, _GIZALIB), Cvoid,
          (Cint, Ptr{Cfloat}, Ptr{Cfloat}), n, xpts, ypts)

function set_line_style(ls::Integer)
    ccall((:giza_set_line_style, _GIZALIB), Cvoid, (Cint,), ls)
end

function get_line_style()
    ls = Ref{Cint}()
    ccall((:giza_get_line_style, _GIZALIB), Cvoid, (Ptr{Cint},), ls)
    return ls[]
end

function set_line_width(lw::Real)
    ccall((:giza_set_line_width, _GIZALIB), Cvoid, (Cdouble,), lw)
end

function get_line_width()
    lw = Ref{Cdouble}()
    ccall((:giza_get_line_width, _GIZALIB), Cvoid, (Ptr{Cdouble},), lw)
    return lw[]
end

function set_line_cap(lc::Integer)
    ccall((:giza_set_line_cap, _GIZALIB), Cvoid, (Cint,), lc)
end

function get_line_cap()
    lc = Ref{Cint}()
    ccall((:giza_get_line_cap, _GIZALIB), Cvoid, (Ptr{Cint},), lc)
    return lc[]
end

for f in (:mark_line, :mark_line_ordered)
    @eval begin
        function $f(maxpts::Integer)
            @assert maxpts ≥ 1
            xpts = Array{Cdouble}(undef, maxpts)
            ypts = Array{Cdouble}(undef, maxpts)
            npts = Ref{Cint}(0)
            ccall(($(string("giza_",f)), _GIZALIB), Cvoid,
                  (Cint, Ptr{Cint}, Ptr{Cdouble}, Ptr{Cdouble}),
                  maxpts, npts, xpts, ypts)
            if npts[] < maxpts
                resize!(xpts, npts[])
                resize!(ypts, npts[])
            end
            return (xpts, ypts)
        end

        function $f(maxpts::Integer, xanc::Real, yanc::Real)
            @assert maxpts ≥ 1
            xpts = Array{Cdouble}(undef, maxpts)
            ypts = Array{Cdouble}(undef, maxpts)
            npts = Ref{Cint}(1)
            xpts[1] = xanc
            ypts[1] = yanc
            ccall(($(string("giza_",f)), _GIZALIB), Cvoid,
                  (Cint, Ptr{Cint}, Ptr{Cdouble}, Ptr{Cdouble}),
                  maxpts, npts, xpts, ypts)
            if npts[] < maxpts
                resize!(xpts, npts[])
                resize!(ypts, npts[])
            end
            return (xpts, ypts)
        end
    end
end

for f in (:mark_points, :mark_points_ordered)
    @eval begin
        function $f(maxpts::Integer, sym::Integer)
            @assert maxpts ≥ 1
            xpts = Array{Cdouble}(undef, maxpts)
            ypts = Array{Cdouble}(undef, maxpts)
            npts = Ref{Cint}(0)
            ccall(($(string("giza_",f)), _GIZALIB), Cvoid,
                  (Cint, Ptr{Cint}, Ptr{Cdouble}, Ptr{Cdouble}, Cint),
                  maxpts, npts, xpts, ypts, sym)
            if npts[] < maxpts
                resize!(xpts, npts[])
                resize!(ypts, npts[])
            end
            return (xpts, ypts)
        end

        function $f(maxpts::Integer, sym::Integer, xanc::Real, yanc::Real)
            @assert maxpts ≥ 1
            xpts = Array{Cdouble}(undef, maxpts)
            ypts = Array{Cdouble}(undef, maxpts)
            npts = Ref{Cint}(1)
            xpts[1] = xanc
            ypts[1] = yanc
            ccall(($(string("giza_",f)), _GIZALIB), Cvoid,
                  (Cint, Ptr{Cint}, Ptr{Cdouble}, Ptr{Cdouble}, Cint),
                  maxpts, npts, xpts, ypts, sym)
            if npts[] < maxpts
                resize!(xpts, npts[])
                resize!(ypts, npts[])
            end
            return (xpts, ypts)
        end
    end
end

function move(x::Real, y::Real)
    ccall((:giza_move, _GIZALIB), Cvoid,
          (Cdouble, Cdouble), x, y)
end

function set_paper_size(units::Integer, width::Real, height::Real)
    ccall((:giza_set_paper_size, _GIZALIB), Cvoid,
          (Cint, Cdouble, Cdouble),
          units, width, height)
end

function get_paper_size(units::Integer)
    width  = Ref{Cdouble}()
    height = Ref{Cdouble}()
    ccall((:giza_get_paper_size, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}),
          units, width, height)
    return (width[], height[])
end

function points(xpts::AbstractVector{Tx},
                ypts::AbstractVector{Ty},
                sym::Integer) where {Tx<:Real, Ty<:Real}
    @assert length(xpts) == length(ypts)
    T = _float(Tx, Ty)
    n = length(xpts)
    _points(n, _array(T, xpts), _array(T, ypts), sym)
end

function _points(n::Integer,
                 xpts::DenseVector{Cdouble},
                 ypts::DenseVector{Cdouble},
                 sym::Integer)
    ccall((:giza_points, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}, Cint),
          n, xpts, ypts, sym)
end

function _points(n::Integer,
                 xpts::DenseVector{Cfloat},
                 ypts::DenseVector{Cfloat},
                 sym::Integer)
    ccall((:giza_points_float, _GIZALIB), Cvoid,
          (Cint, Ptr{Cfloat}, Ptr{Cfloat}, Cint),
          n, xpts, ypts, sym)
end

function single_point(x::Real, y::Real, sym::Integer)
    ccall((:giza_single_point, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cint),
          x, y, sym)
end

function polygon(xpts::AbstractVector{Tx},
                 ypts::AbstractVector{Ty}) where {Tx<:Real, Ty<:Real}
    @assert length(xpts) == length(ypts)
    T = _float(Tx, Ty)
    n = length(xpts)
    _polygon(n, _array(T, xpts), _array(T, ypts))
end

function _polygon(n::Integer,
                  xpts::DenseVector{Cdouble},
                  ypts::DenseVector{Cdouble})
    ccall((:giza_polygon, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}),
          n, xpts, ypts)
end

function _polygon(n::Integer,
                  xpts::DenseVector{Cfloat},
                  ypts::DenseVector{Cfloat})
    ccall((:giza_polygon_float, _GIZALIB), Cvoid,
          (Cint, Ptr{Cfloat}, Ptr{Cfloat}),
          n, xpts, ypts)
end

function text(x::Real, y::Real, text::AbstractString)
    ccall((:giza_text, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cstring),
          x, y, text)
end

function ptext(x::Real, y::Real, angle::Real, just::Real, text::AbstractString)
    ccall((:giza_ptext, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble, Cstring),
          x, y, angle, just, text)
end

function qtext(x::Real, y::Real, angle::Real, just::Real, text::AbstractString)
    xbox = Ref{NTuple{4,Cdouble}}()
    ybox = Ref{NTuple{4,Cdouble}}()
    ccall((:giza_qtext, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble, Cstring,
           Ptr{Cdouble}, Ptr{Cdouble}),
          x, y, angle, just, text, xbox, ybox)
    return (xbox[], ybox[])
end

function qtextlen(units::Integer, text::AbstractString)
    xlen = Ref{Cdouble}()
    ylen = Ref{Cdouble}()
    ccall((:giza_qtextlen, _GIZALIB), Cvoid,
          (Cint, Cstring, Ptr{Cdouble}, Ptr{Cdouble}),
          units, text, xlen, ylen)
    return (xlen[], ylen[])
end

function rectangle(x1::Real, x2::Real, y1::Real, y2::Real)
    ccall((:giza_rectangle, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble),
          x1, x2, y1, y2)
end

function rectangle_rounded(x1::Real, x2::Real, y1::Real, y2::Real, r::Real)
    ccall((:giza_rectangle_rounded, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble, Cdouble),
          x1, x2, y1, y2, r)
end

for f in (:begin_autolog,
          :begin_buffer,
          :change_page,
          :close_device,
          :draw_background,
          :end_autolog,
          :end_buffer,
          :flush_buffer,
          :flush_device,
          :print_device_list,
          :print_id,
          :restore,
          :restore_colour_table,
          :save,
          :save_colour_table,
          :set_viewport_default,
          :start_prompting,
          :start_warnings,
          :stop_prompting,
          :stop_warnings)
    @eval $f() =
        ccall(($(string("giza_",f)), _GIZALIB), Cvoid, ())
end

for f in (:render, :render_transparent, :render_gray)
    _f = Symbol(:_,f)
    @eval begin
        function $f(arr::AbstractMatrix{Tarr},
                    vmin::Real,
                    vmax::Real,
                    extend::Integer,
                    aff::AbstractArray{Taff,N}) where {Tarr<:Real,
                                                       Taff<:Real, N}
            @assert length(aff) == 6
            T = _float(Tarr, Taff)
            xdim, ydim = size(arr)
            $_f(xdim, ydim, _array(T, arr), 0, xdim - 1, 0, ydim - 1,
                convert(T, vmin), convert(T, vmax), extend,
                _array(T, aff))
        end

        function $_f(xdim::Integer,
                     ydim::Integer,
                     arr::DenseMatrix{Cdouble},
                     i1::Integer,
                     i2::Integer,
                     j1::Integer,
                     j2::Integer,
                     vmin::Cdouble,
                     vmax::Cdouble,
                     extend::Integer,
                     aff::DenseArray{Cdouble,N}) where {N}
            ccall(($(string("giza_",f)), _GIZALIB), Cvoid,
                  (Cint, Cint, Ptr{Cdouble}, Cint, Cint, Cint, Cint,
                   Cdouble, Cdouble, Cint, Ptr{Cdouble}),
                  xdim, ydim, arr, i1, i2, j1, j2, vmin, vmax, extend, aff)
        end

        function $_f(xdim::Integer,
                     ydim::Integer,
                     arr::DenseMatrix{Cfloat},
                     i1::Integer,
                     i2::Integer,
                     j1::Integer,
                     j2::Integer,
                     vmin::Cfloat,
                     vmax::Cfloat,
                     extend::Integer,
                     aff::DenseArray{Cfloat,N}) where {N}
            ccall(($(string("giza_",f,"_float")), _GIZALIB), Cvoid,
                  (Cint, Cint, Ptr{Cfloat}, Cint, Cint, Cint, Cint,
                   Cfloat, Cfloat, Cint, Ptr{Cfloat}),
                  xdim, ydim, arr, i1, i2, j1, j2, vmin, vmax, extend, aff)
        end

    end
end

function render_alpha(arr::AbstractMatrix{Tarr},
                      alpha::AbstractMatrix{Talpha},
                      vmin::Real,
                      vmax::Real,
                      extend::Integer,
                      aff::AbstractArray{Taff,N}) where {Tarr<:Real,
                                                         Talpha<:Real,
                                                         Taff<:Real, N}
    @assert size(arr) == size(alpha)
    @assert length(aff) == 6
    T = _float(Tarr, Talpha, Taff)
    xdim, ydim = size(arr)
    _render_alpha(xdim, ydim, _array(T, arr),  _array(T, alpha),
                  0, xdim - 1, 0, ydim - 1,
                  convert(T, vmin), convert(T, vmax),
                  extend, _array(T, aff))
end

function _render_alpha(xdim::Integer,
                       ydim::Integer,
                       arr::DenseMatrix{Cdouble},
                       alpha::DenseMatrix{Cdouble},
                       i1::Integer,
                       i2::Integer,
                       j1::Integer,
                       j2::Integer,
                       vmin::Cdouble,
                       vmax::Cdouble,
                       extend::Integer,
                       aff::DenseArray{Cdouble,N}) where {N}
    ccall((:giza_render_alpha, _GIZALIB), Cvoid,
          (Cint, Cint, Ptr{Cdouble}, Ptr{Cdouble}, Cint, Cint, Cint, Cint,
           Cdouble, Cdouble, Cint, Ptr{Cdouble}),
          xdim, ydim, arr, alpha, i1, i2, j1, j2, vmin, vmax, extend, aff)
end

function _render_alpha(xdim::Integer,
                       ydim::Integer,
                       arr::DenseArray{Cfloat},
                       alpha::DenseArray{Cfloat},
                       i1::Integer,
                       i2::Integer,
                       j1::Integer,
                       j2::Integer,
                       vmin::Cfloat,
                       vmax::Cfloat,
                       extend::Integer,
                       aff::DenseArray{Cfloat,N}) where {N}
    ccall((:giza_render_alpha_float, _GIZALIB), Cvoid,
          (Cint, Cint, Ptr{Cfloat}, Ptr{Cfloat}, Cint, Cint, Cint, Cint,
           Cfloat, Cfloat, Cint, Ptr{Cfloat}),
          xdim, ydim, arr, alpha, i1, i2, j1, j2, vmin, vmax, extend, aff)
end

function draw_pixels(arr::AbstractMatrix{<:Integer},
                     xmin::Real,
                     xmax::Real,
                     ymin::Real,
                     ymax::Real,
                     extend::Integer)
    xdim, ydim = size(arr)
    ccall((:giza_draw_pixels, _GIZALIB), Cvoid,
          (Cint, Cint, Ptr{Cint}, Cint, Cint, Cint, Cint,
           Cdouble, Cdouble, Cdouble, Cdouble, Cint),
          xdim, ydim, _array(Cint, arr), 0, xdim - 1, 0, ydim - 1,
          xmin, xmax, ymin, ymax, extend)
end

function nice_number(x::Real)
    nsub = Ref{Cint}()
    xrnd = ccall((:giza_round, _GIZALIB), Cdouble,
                 (Cdouble, Ptr{Cint}), x, nsub)
    return (xrnd, nsub[])
end

function subpanel(nx::Integer, ny::Integer)
    ccall((:giza_subpanel, _GIZALIB), Cvoid,
          (Cint, Cint), nx, ny)
end

function set_panel(ix::Integer, iy::Integer)
    ccall((:giza_set_panel, _GIZALIB), Cvoid,
          (Cint, Cint), ix, iy)
end

function get_panel()
    ix = Ref{Cint}()
    iy = Ref{Cint}()
    ccall((:giza_get_panel, _GIZALIB), Cvoid,
          (Ptr{Cint}, Ptr{Cint}), ix, iy)
    return (ix[], iy[])
end

for f in (:set_font,
          :set_font_bold,
          :set_font_bold_italic,
          :set_font_italic)
    @eval function $f(font::AbstractString)
        ccall(($(string("giza_",f)), _GIZALIB), Cvoid,
              (Cstring,), font)
    end
end

function get_font()
    len = 200
    buf = Array{UInt8}(undef, len)
    ccall((:giza_get_font, _GIZALIB), Cvoid,
          (Ptr{UInt8}, Cint), buf, len)
    buf[len] = 0
    return unsafe_string(pointer(buf))
end

function set_text_background(ci::Integer)
    ccall((:giza_set_text_background, _GIZALIB), Cvoid, (Cint,), ci)
end

function get_text_background()
    ci = Ref{Cint}()
    ccall((:giza_get_text_background, _GIZALIB), Cvoid, (Ptr{Cint},), ci)
    return ci[]
end

function vector(vx::AbstractMatrix{Tx},
                vy::AbstractMatrix{Ty},
                scale::Real,
                position::Integer,
                aff::AbstractArray{Ta,N},
                blank::Real) where {Tx<:Real, Ty<:Real, Ta<:Real, N}
    @assert size(vx) == size(vy)
    @assert length(aff) == 6
    T = _float(Tx, Ty, Ta)
    xdim, ydim = size(vx)
    _vector(xdim, ydim, _array(T, vx), _array(T, vy),
            0, xdim - 1, 0, ydim - 1,
            scale, position, aff, blank)
end

function _vector(xdim::Integer,
                 ydim::Integer,
                 vx::DenseVector{Cdouble},
                 vy::DenseVector{Cdouble},
                 i1::Integer,
                 i2::Integer,
                 j1::Integer,
                 j2::Integer,
                 scale::Real,
                 position::Integer, # FIXME:
                 aff::DenseArray{Cdouble,N},
                 blank::Real) where {N}
    ccall((:giza_vector, _GIZALIB), Cvoid,
          (Cint, Cint, Ptr{Cdouble}, Ptr{Cdouble},
           Cint, Cint, Cint, Cint, Cdouble, Cint,
           Ptr{Cdouble}, Cdouble),
          xdim, ydim, _array(T, vx), _array(T, vy),
          0, xdim - 1, 0, ydim - 1, scale, position, aff, blank)
end

function _vector(xdim::Integer,
                 ydim::Integer,
                 vx::DenseVector{Cfloat},
                 vy::DenseVector{Cfloat},
                 i1::Integer,
                 i2::Integer,
                 j1::Integer,
                 j2::Integer,
                 scale::Real,
                 position::Integer, # FIXME:
                 aff::DenseArray{Cfloat,N},
                 blank::Real) where {N}
    ccall((:giza_vector_float, _GIZALIB), Cvoid,
          (Cint, Cint, Ptr{Cfloat}, Ptr{Cfloat},
           Cint, Cint, Cint, Cint, Cfloat, Cint,
           Ptr{Cfloat}, Cfloat),
          xdim, ydim, _array(T, vx), _array(T, vy),
          0, xdim - 1, 0, ydim - 1, scale, position, aff, blank)
end

function set_viewport(left::Real, right::Real, bottom::Real, top::Real)
    ccall((:giza_set_viewport, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble),
          left, right, bottom, top)
end

function get_viewport(units::Integer)
    left   = Ref{Cdouble}()
    right  = Ref{Cdouble}()
    bottom = Ref{Cdouble}()
    top    = Ref{Cdouble}()
    ccall((:giza_get_viewport, _GIZALIB), Cvoid,
          (Cint, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}),
          units, left, right, bottom, top)
    return (left[], right[], bottom[], top[])
end

function set_viewport_inches(left::Real, right::Real, bottom::Real, top::Real)
    ccall((:giza_set_viewport_inches, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble),
          left, right, bottom, top)
end

function version()
    major = Ref{Cint}()
    minor = Ref{Cint}()
    micro = Ref{Cint}()
    ccall((:giza_version, _GIZALIB), Cvoid,
          (Ptr{Cint}, Ptr{Cint}, Ptr{Cint}),
          major, minor, micro)
    return (major[], minor[], micro[])
end

function set_window(x1::Real, x2::Real, y1::Real, y2::Real)
    ccall((:giza_set_window, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble),
          x1, x2, y1, y2)
end

function set_window_equal_scale(x1::Real, x2::Real, y1::Real, y2::Real)
    ccall((:giza_set_window_equal_scale, _GIZALIB), Cvoid,
          (Cdouble, Cdouble, Cdouble, Cdouble),
          x1, x2, y1, y2)
end

function get_window()
    x1 = Ref{Cdouble}()
    x2 = Ref{Cdouble}()
    y1 = Ref{Cdouble}()
    y2 = Ref{Cdouble}()
    ccall((:giza_get_window, _GIZALIB), Cvoid,
          (Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}),
          x1, x2, y1, y2)
    return (x1[], x2[], y1[], y2[])
end

end # module


text:
# 2bit Writer
# ===========

"""
    TwoBitWriter(output::IO, names::AbstractVector)

Create a data writer of the 2bit file format.

# Arguments
* `output`: data sink
* `names`: a vector of sequence names written to `output`
"""
type TwoBitWriter{T<:IO} <: Bio.IO.AbstractWriter
    # output stream
    output::T

    # sequence names
    names::Vector{String}

    # bit vector to check if each sequence is already written or not
    written::BitVector
end

function TwoBitWriter(output::IO, names::AbstractVector)
    writer = TwoBitWriter(output, names, falses(length(names)))
    write_header(writer)
    write_index(writer)
    return writer
end

function Bio.IO.stream(writer::TwoBitWriter)
    return writer.output
end

function Base.close(writer::TwoBitWriter)
    if !all(writer.written)
        error("one or more sequences are not written")
    end
    close(writer.output)
end

function write_header(writer::TwoBitWriter)
    output = writer.output
    n = 0
    n += write(output, 0x1A412743)
    n += write(output, UInt32(0))
    n += write(output, UInt32(length(writer.names)))
    n += write(output, UInt32(0))
    return n
end

function write_index(writer::TwoBitWriter)
    output = writer.output
    n = 0
    for name in writer.names
        n += write(output, UInt8(length(name)))
        n += write(output, name)
        n += write(output, UInt32(0))
    end
    return n
end

# Update the file offset of a sequence in the index section.
function update_offset(writer::TwoBitWriter, seqname, seqoffset)
    @assert seqname ∈ writer.names
    output = writer.output
    offset = 16
    for name in writer.names
        offset += sizeof(UInt8) + length(name)
        if name == seqname
            old = position(output)
            seek(output, offset)
            write(output, UInt32(seqoffset))
            seek(output, old)
            return
        end
        offset += sizeof(UInt32)
    end
end

function Base.write(writer::TwoBitWriter, record::SeqRecord)
    i = findfirst(writer.names, record.name)
    if i == 0
        error("sequence \"", record.name, "\" doesn't exist in the writing list")
    elseif writer.written[i]
        error("sequence \"", record.name, "\" is already written")
    end

    output = writer.output
    update_offset(writer, record.name, position(output))

    n = 0
    n += write(output, UInt32(length(record.seq)))
    n += write_n_blocks(output, record.seq)
    n += write_masked_blocks(output, record.metadata)
    n += write(output, UInt32(0))  # reserved bytes
    n += write_twobit_sequence(output, record.seq)

    writer.written[i] = true
    return n
end

function make_n_blocks(seq)
    starts = UInt32[]
    sizes = UInt32[]
    i = 1
    while i ≤ endof(seq)
        nt = seq[i]
        if nt == DNA_N
            start = i - 1  # 0-based index
            push!(starts, start)
            while i ≤ endof(seq) && seq[i] == DNA_N
                i += 1
            end
            push!(sizes, (i - 1) - start)
        elseif isambiguous(nt)
            error("ambiguous nucleotide except N is not supported")
        else
            i += 1
        end
    end
    return starts, sizes
end

function write_n_blocks(output, seq)
    blockstarts, blocksizes = make_n_blocks(seq)
    @assert length(blockstarts) == length(blocksizes)
    n = 0
    n += write(output, UInt32(length(blockstarts)))
    n += write(output, blockstarts)
    n += write(output, blocksizes)
    return n
end

function write_masked_blocks(output, metadata)
    n = 0
    if isa(metadata, Vector{UnitRange{Int}})
        n += write(output, UInt32(length(metadata)))
        for mblock in metadata
            n += write(output, UInt32(first(mblock) - 1))  # 0-based
        end
        for mblock in metadata
            n += write(output, UInt32(length(mblock)))
        end
    elseif metadata === nothing
        n += write(output, UInt32(0))
    else
        error("metadata is not serializable in the 2bit file format")
    end
    return n
end

function write_twobit_sequence(output, seq)
    n = 0
    i = 4
    while i ≤ endof(seq)
        x::UInt8 = 0
        x |= nuc2twobit(seq[i-3]) << 6
        x |= nuc2twobit(seq[i-2]) << 4
        x |= nuc2twobit(seq[i-1]) << 2
        x |= nuc2twobit(seq[i-0]) << 0
        n += write(output, x)
        i += 4
    end
    r = length(seq) % 4
    if r > 0
        let x::UInt8 = 0
            i = endof(seq) - r + 1
            while i ≤ endof(seq)
                x = x << 2 | nuc2twobit(seq[i])
                i += 1
            end
            x <<= (4 - r) * 2
            n += write(output, x)
        end
    end
    return n
end

function nuc2twobit(nt::DNA)
    return (
        nt == DNA_A ? 0b10 :
        nt == DNA_C ? 0b01 :
        nt == DNA_G ? 0b11 :
        nt == DNA_T ? 0b00 :
        nt == DNA_N ? 0b00 : error())
end


text:
module TestRLEDF

using Test
using RLEVectors

@testset begin

    # Initialization
    x = RLEDataFrame([RLEVector([1, 1, 2]), RLEVector([2, 2, 2])], [:a, :b])
    @test isa(x, RLEDataFrame)
    @test length(x) == 2
    @test names(x) == [:a, :b]
    io = IOBuffer()
    @test typeof(show(io, x)) == Nothing # At least test that show does not give error
    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    @test isa(z, RLEDataFrame)
    @test length(z) == 3
    @test names(z) == [:a, :b, :c]
    @test columns(z) == [RLEVector([5, 2, 2]), RLEVector([4, 4, 4]), RLEVector([3, 2, 1])]
    @test names(z) == [:a, :b, :c]
    @test_throws ArgumentError RLEDataFrame([RLEVector([1])], [:a, :b])
    @test_throws ArgumentError RLEDataFrame([RLEVector([1]), RLEVector([2, 3])], [:a, :b])

    # Getting and setting
    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    @test z[:] == z
    @test z[2] == RLEVector([4, 4, 4])
    @test z[:c] == RLEVector([3, 2, 1])
    @test z[2:3] == RLEDataFrame(b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    @test z[[:b, :c]] == RLEDataFrame(b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    @test z[[2], [:b, :c]] == RLEDataFrame(b = RLEVector([4]), c = RLEVector([2]))
    @test z[2, [:b, :c]] == RLEDataFrame(b = RLEVector([4]), c = RLEVector([2]))
    @test z[2, 3] == 2
    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    @test z[1, 3] == 3

    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    z[2] = RLEVector([9, 8, 9])
    @test z ==
          RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([9, 8, 9]), c = RLEVector([3, 2, 1]))
    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    z[:b] = RLEVector([7, 2, 1])
    @test z ==
          RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([7, 2, 1]), c = RLEVector([3, 2, 1]))
    z[:b] = RLEVector([7, 2, 1])
    @test z ==
          RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([7, 2, 1]), c = RLEVector([3, 2, 1]))
    z[:d] = RLEVector([2, 1, 2])
    @test z == RLEDataFrame(
        a = RLEVector([5, 2, 2]),
        b = RLEVector([7, 2, 1]),
        c = RLEVector([3, 2, 1]),
        d = RLEVector([2, 1, 2]),
    )
    @test_throws ArgumentError z[1] = RLEVector([2])
    @test_throws BoundsError z[14] = RLEVector([2, 2, 1])
    @test_throws ArgumentError z[:a] = RLEVector([2, 4, 5, 4, 3, 4])

    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    z[1:2, :b] = [5, 5]
    @test z ==
          RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([5, 5, 4]), c = RLEVector([3, 2, 1]))

    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    z[1:2, [2]] = [5, 5]
    @test z ==
          RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([5, 5, 4]), c = RLEVector([3, 2, 1]))

    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    z[1, 3] = 12
    @test z == RLEDataFrame(
        a = RLEVector([5, 2, 2]),
        b = RLEVector([4, 4, 4]),
        c = RLEVector([12, 2, 1]),
    )

    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    z[1, [2, 3]] = 12
    @test z == RLEDataFrame(
        a = RLEVector([5, 2, 2]),
        b = RLEVector([12, 4, 4]),
        c = RLEVector([12, 2, 1]),
    )

    z = RLEDataFrame(a = RLEVector([5, 2, 2]), b = RLEVector([4, 4, 4]), c = RLEVector([3, 2, 1]))
    #@test rowMeans(z) == [7,8/3,7/3]
    #@test rowSums(z) == [12,8,7]
    #@test rowMedians(z) == [5,4,4]
    @test colMeans(z) == [3, 4, 2]
    @test colSums(z) == [9, 12, 6]
    @test colMedians(z) == [2.0, 4.0, 2.0]

    zm = convert(Matrix, z)
    @test rowMeans(zm) == [12 / 3, 8 / 3, 7 / 3]
    @test rowSums(zm) == [12, 8, 7]
    @test rowMedians(zm) == [4, 2, 2]
    @test colMeans(zm) == [3, 4, 2]
    @test colSums(zm) == [9, 12, 6]
    @test colMedians(zm) == [2.0, 4.0, 2.0]

end # testset
end # module


text:
using NiLang, NiLang.AD
using CuArrays, CUDAnative, GPUArrays
using BenchmarkTools

@i @inline function ⊖(CUDAnative.pow)(out!::GVar{T}, x::GVar{T}, n::GVar) where T
    ⊖(CUDAnative.pow)(value(out!), value(x), value(n))

    # grad x
    @routine @invcheckoff begin
        @zeros T anc1 anc2 anc3 jac1 jac2

        DEC(value(n))
        anc1 += CUDAnative.pow(value(x), value(n))
        INC(value(n))
        jac1 += anc1 * value(n)

        # get grad of n
        anc2 += log(value(x))
        anc3 += CUDAnative.pow(value(x), value(n))
        jac2 += anc3*anc2
    end
    grad(x) += grad(out!) * jac1
    grad(n) += grad(out!) * jac2
    ~@routine
end

@i @inline function ⊖(CUDAnative.pow)(out!::GVar{T}, x::GVar, n) where T
    ⊖(CUDAnative.pow)(value(out!), value(x), n)
    @routine @invcheckoff begin
        anc1 ← zero(value(x))
        jac ← zero(value(x))

        DEC(value(n))
        anc1 += CUDAnative.pow(value(x), n)
        INC(value(n))
        jac += anc1 * n
    end
    grad(x) += grad(out!) * jac
    ~@routine
end

@i @inline function ⊖(CUDAnative.pow)(out!::GVar{T}, x, n::GVar) where T
    ⊖(CUDAnative.pow)(value(out!), x, value(n))
    # get jac of n
    @routine @invcheckoff begin
        anc1 ← zero(x)
        anc2 ← zero(x)
        jac ← zero(x)

        anc1 += log(x)
        anc2 += CUDAnative.pow(x, value(n))
        jac += anc1*anc2
    end
    grad(n) += grad(out!) * jac
    ~@routine
end

# You need to replace all "^" operations in `ibessel` with `CUDAnative.pow`.
# Please remember to turn invertiblity check off, because error handling is not supported in a cuda thread.
# Function `i_dirtymul` and `i_factorial` are not changed.

@i function ibesselj(out!, ν, z; atol=1e-8)
    @routine @invcheckoff begin
        k ← 0
        fact_nu ← zero(ν)
        halfz ← zero(z)
        halfz_power_nu ← zero(z)
        halfz_power_2 ← zero(z)
        out_anc ← zero(z)
        anc1 ← zero(z)
        anc2 ← zero(z)
        anc3 ← zero(z)
        anc4 ← zero(z)
        anc5 ← zero(z)

        halfz += z / 2
        halfz_power_nu += CUDAnative.pow(halfz, ν)
        halfz_power_2 += CUDAnative.pow(halfz, 2)
        i_factorial(fact_nu, ν)

        anc1 += halfz_power_nu/fact_nu
        out_anc += anc1
        @from k==0 while abs(unwrap(anc1)) > atol && abs(unwrap(anc4)) < atol
            INC(k)
            @routine begin
                anc5 += k
                anc5 += ν
                anc2 -= k * anc5
                anc3 += halfz_power_2 / anc2
            end
            i_dirtymul(anc1, anc3, anc4)
            out_anc += anc1
            ~@routine
        end
    end
    out! += out_anc
    ~@routine
end

# Define your reversible kernel function that calls the reversible bessel function

@i function ibesselj_kernel(out!, ν, z, atol)
    i ← (blockIdx().x-1) * blockDim().x + threadIdx().x
    @inbounds ibesselj(out![i], ν, z[i]; atol=atol)
    @invcheckoff i → (blockIdx().x-1) * blockDim().x + threadIdx().x
end

# To launch this reversible kernel, you also need a reversible host function.

@i function ibesselj(out!::CuVector, ν, z::CuVector; atol=1e-8)
   XY ← GPUArrays.thread_blocks_heuristic(length(out!))
   @cuda threads=XY.:1 blocks=XY.:2 ibesselj_kernel(out!, ν, z, atol)
   @invcheckoff XY → GPUArrays.thread_blocks_heuristic(length(out!))
end

# To test this function, we first define input parameters `a` and output `out!`
N = 4096
T = Float64
a = CuArray(ones(T, N))
out! = CuArray(zeros(T, N))

# We wrap the output with a randomly initialized gradient field, suppose we get the gradients from a virtual loss function.
# Also, we need to initialize an empty gradient field for elements in input cuda tensor `a`.
out! = ibesselj(out!, 2, GVar.(a))[1]
out_g! = GVar.(out!, CuArray(ones(T, N)))
a_g = GVar.(a)

# Call the inverse program, the multiple dispatch will drive you to the goal.
println("Benchmarking NiLang on CUDA, N = $N, T = $T")
display(@benchmark CuArrays.@sync (~ibesselj)($out_g!, 2, $a_g))


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

"""
    delete_record(event_time, feature_group_name, record_identifier_value_as_string)
    delete_record(event_time, feature_group_name, record_identifier_value_as_string, params::Dict{String,<:Any})

Deletes a Record from a FeatureGroup. A new record will show up in the OfflineStore when
the DeleteRecord API is called. This record will have a value of True in the is_deleted
column.

# Arguments
- `event_time`: Timestamp indicating when the deletion event occurred. EventTime can be
  used to query data at a certain point in time.
- `feature_group_name`: The name of the feature group to delete the record from.
- `record_identifier_value_as_string`: The value for the RecordIdentifier that uniquely
  identifies the record, in string format.

"""
delete_record(EventTime, FeatureGroupName, RecordIdentifierValueAsString; aws_config::AbstractAWSConfig=global_aws_config()) = sagemaker_featurestore_runtime("DELETE", "/FeatureGroup/$(FeatureGroupName)", Dict{String, Any}("EventTime"=>EventTime, "RecordIdentifierValueAsString"=>RecordIdentifierValueAsString); aws_config=aws_config)
delete_record(EventTime, FeatureGroupName, RecordIdentifierValueAsString, params::AbstractDict{String, <:Any}; aws_config::AbstractAWSConfig=global_aws_config()) = sagemaker_featurestore_runtime("DELETE", "/FeatureGroup/$(FeatureGroupName)", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("EventTime"=>EventTime, "RecordIdentifierValueAsString"=>RecordIdentifierValueAsString), params)); aws_config=aws_config)

"""
    get_record(feature_group_name, record_identifier_value_as_string)
    get_record(feature_group_name, record_identifier_value_as_string, params::Dict{String,<:Any})

Use for OnlineStore serving from a FeatureStore. Only the latest records stored in the
OnlineStore can be retrieved. If no Record with RecordIdentifierValue is found, then an
empty result is returned.

# Arguments
- `feature_group_name`: The name of the feature group in which you want to put the records.
- `record_identifier_value_as_string`: The value that corresponds to RecordIdentifier type
  and uniquely identifies the record in the FeatureGroup.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"FeatureName"`: List of names of Features to be retrieved. If not specified, the latest
  value for all the Features are returned.
"""
get_record(FeatureGroupName, RecordIdentifierValueAsString; aws_config::AbstractAWSConfig=global_aws_config()) = sagemaker_featurestore_runtime("GET", "/FeatureGroup/$(FeatureGroupName)", Dict{String, Any}("RecordIdentifierValueAsString"=>RecordIdentifierValueAsString); aws_config=aws_config)
get_record(FeatureGroupName, RecordIdentifierValueAsString, params::AbstractDict{String, <:Any}; aws_config::AbstractAWSConfig=global_aws_config()) = sagemaker_featurestore_runtime("GET", "/FeatureGroup/$(FeatureGroupName)", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("RecordIdentifierValueAsString"=>RecordIdentifierValueAsString), params)); aws_config=aws_config)

"""
    put_record(feature_group_name, record)
    put_record(feature_group_name, record, params::Dict{String,<:Any})

Used for data ingestion into the FeatureStore. The PutRecord API writes to both the
OnlineStore and OfflineStore. If the record is the latest record for the recordIdentifier,
the record is written to both the OnlineStore and OfflineStore. If the record is a historic
record, it is written only to the OfflineStore.

# Arguments
- `feature_group_name`: The name of the feature group that you want to insert the record
  into.
- `record`: List of FeatureValues to be inserted. This will be a full over-write. If you
  only want to update few of the feature values, do the following:   Use GetRecord to
  retrieve the latest record.   Update the record returned from GetRecord.    Use PutRecord
  to update feature values.

"""
put_record(FeatureGroupName, Record; aws_config::AbstractAWSConfig=global_aws_config()) = sagemaker_featurestore_runtime("PUT", "/FeatureGroup/$(FeatureGroupName)", Dict{String, Any}("Record"=>Record); aws_config=aws_config)
put_record(FeatureGroupName, Record, params::AbstractDict{String, <:Any}; aws_config::AbstractAWSConfig=global_aws_config()) = sagemaker_featurestore_runtime("PUT", "/FeatureGroup/$(FeatureGroupName)", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("Record"=>Record), params)); aws_config=aws_config)


text:
# Node tests
a = Node("a",Intervals([(0.0,1.0),(1.0,2.0)]))
b = Node("a",Intervals([(0.0,1.0),(1.0,2.0)]))
c = Node("c",Intervals([(0.000001,0.99)]))

@test a==b
@test b==a
@test a==a
@test a!=c
@test b!=c
@test c ⊆ a
@test a ⊆ b
@test a ⊈ c
@test b ⊈ c
@test a ⊆ a
@test (0.5,0.99) ⊆ a
@test (0.6,1.5) ⊆ b
@test (-0.001,0.9) ⊈ c
@test (0.0,0.99) ⊆ c
@test 0.5 ∈ a
@test 0.98 ∈ c
@test 2.5 ∉ a
@test 2.0 ∉ c
@test (a ∩ c).list == [(0.000001,0.99)]
@test (a ∪ c).list == [(0.0,2.0)]

a = Node("a",Intervals([(0.0,1.0),(1.0,2.0)]))
b = Node("b",Intervals([(0.0,1.5),(1.6,2.5)]))
c = Node("c",Intervals([(1.0,1.6)]))
d = Node("d",Intervals([(0.3,1.2)]))
L = Vector([a,b,c])
@test a ∈ L
@test b ∈ L
@test d ∉ L
@test c ∈ L
@test get_idx(a,L) == [1]
@test get_idx(b,L) == [2]
@test get_idx(d,L) == []
@test get_idx("a",L) == [1]
@test get_idx("b",L) == [2]
@test get_idx("d",L) == []
@test get_idx(1.0,L) == [1,2,3]
@test get_idx(1.61,L) == [1,2]
@test get_idx(2.4,L) == [2]
@test get_idx(2.6,L) == []
L2 = Vector([a,b])
@test L2 ⊆ L
@test L ∩ L2 == [a,b]
@test L ∪ L2 == [a,b,c]
L2 = Vector([a,b,d])
@test L2 ⊈ L
@test L ∩ L2 == [a,b]
@test L ∪ L2 == [a,b,c,d]
@test L2 ∪ L == [a,b,d,c]
L2 = Vector([])
@test L2 ⊆ L
@test L ∩ L2 == []
@test L ∪ L2 == L
L2 = [Node("a",Intervals([(0.3,0.9)])),
      Node("x",Intervals([(1.3,1.5)])),
      Node("z",Intervals([(1.2,1.3)]))]
@test all(L2 .⊆ L)
L2 = [Node("a",Intervals([(0.3,0.9)])),
      Node("b",Intervals([(1.3,1.9)])),
      Node("z",Intervals([(1.2,1.3)]))]
@test !all(L2 .⊆ L)

a1 = Node("a",Intervals([(0.5,1.0),(1.5,2.5)]))
a2 = Node("a",Intervals([(0.0,0.45),(0.6,1.0),(1.1,1.4),(2.3,3.0)]))
a3 = merge(a1,a2)
@test a3.name == "a"
@test a3.presence == Intervals([(0.0,0.45),(0.5,1.0),(1.1,1.4),(1.5,3.0)])
merge!(a1,a2)
@test a1.name == "a"
@test a1.presence == Intervals([(0.0,0.45),(0.5,1.0),(1.1,1.4),(1.5,3.0)])