text:
@enum AggOperations begin

    MEAN = 0
    MEDIAN = 1
    EACH = 2

end

AggType = Union{AggOperations, Integer}

@enum FilterTypes begin

    BUTTERWORTH = 0
    CHEBYSHEV_TYPE_1 = 1
    BESSEL = 2

end

FiltType = Union{FilterTypes, Integer}

@enum WindowOperations begin

    NO_WINDOW = 0
    HANNING = 1
    HAMMING = 2
    BLACKMAN_HARRIS = 3

end

WinType = Union{WindowOperations, Integer}

@enum DetrendOperations begin

    NO_DETREND = 0
    CONSTANT = 1
    LINEAR = 2

end

DetType = Union{DetrendOperations, Integer}

@enum NoiseTypes begin

    FIFTY = 0
    SIXTY = 1
    FIFTY_AND_SIXTY = 2

end

EnvNoiseType = Union{NoiseTypes, Integer}

@enum WaveletDenoisingTypes begin

    VISUSHRINK = 0
    SURESHRINK = 1

end

WaveletDenoisingType = Union{WaveletDenoisingTypes, Integer}

@enum ThresholdTypes begin

    SOFT = 0
    HARD = 1

end

ThresholdType = Union{ThresholdTypes, Integer}

@enum WaveletExtensionTypes begin

    SYMMETRIC = 0
    PERIODIC = 1

end

WaveletExtensionType = Union{WaveletExtensionTypes, Integer}

@enum NoiseEstimationLevelTypes begin

    FIRST_LEVEL = 0
    ALL_LEVELS = 1

end

NoiseEstimationLevelType = Union{NoiseEstimationLevelTypes, Integer}

@enum WaveletTypes begin

    HAAR = 0
    DB1 = 1
    DB2 = 2
    DB3 = 3
    DB4 = 4
    DB5 = 5
    DB6 = 6
    DB7 = 7
    DB8 = 8
    DB9 = 9
    DB10 = 10
    DB11 = 11
    DB12 = 12
    DB13 = 13
    DB14 = 14
    DB15 = 15
    BIOR1_1 = 16
    BIOR1_3 = 17
    BIOR1_5 = 18
    BIOR2_2 = 19
    BIOR2_4 = 20
    BIOR2_6 = 21
    BIOR2_8 = 22
    BIOR3_1 = 23
    BIOR3_3 = 24
    BIOR3_5 = 25
    BIOR3_7 = 26
    BIOR3_9 = 27
    BIOR4_4 = 28
    BIOR5_5 = 29
    BIOR6_8 = 30
    COIF1 = 31
    COIF2 = 32
    COIF3 = 33
    COIF4 = 34
    COIF5 = 35
    SYM2 = 36
    SYM3 = 37
    SYM4 = 38
    SYM5 = 39
    SYM6 = 40
    SYM7 = 41
    SYM8 = 42
    SYM9 = 43
    SYM10 = 44

end

WaveletType = Union{WaveletTypes, Integer}

@brainflow_rethrow function perform_lowpass(data, sampling_rate::Integer, cutoff::Float64, order::Integer,
    filter_type::FiltType, ripple::Float64)
    ccall((:perform_lowpass, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Float64, Cint, Cint, Float64),
            data, length(data), Int32(sampling_rate), Float64(cutoff), Int32(order), Int32(filter_type), Float64(ripple))
    return
end

@brainflow_rethrow function perform_highpass(data, sampling_rate::Integer, cutoff::Float64, order::Integer,
    filter_type::FiltType, ripple::Float64
)
    ccall((:perform_highpass, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Float64, Cint, Cint, Float64),
            data, length(data), Int32(sampling_rate), Float64(cutoff), Int32(order), Int32(filter_type), Float64(ripple))
    return
end

@brainflow_rethrow function perform_bandpass(data, sampling_rate::Integer, start_freq::Float64,
    stop_freq::Float64, order::Integer, filter_type::FiltType, ripple::Float64)
    ccall((:perform_bandpass, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Float64, Float64, Cint, Cint, Float64),
            data, length(data), Int32(sampling_rate), Float64(start_freq), Float64(stop_freq), Int32(order), Int32(filter_type), Float64(ripple))
    return
end

@brainflow_rethrow function perform_bandstop(data, sampling_rate::Integer, start_freq::Float64,
    stop_freq::Float64, order::Integer, filter_type::FiltType, ripple::Float64)
    ccall((:perform_bandstop, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Float64, Float64, Cint, Cint, Float64),
            data, length(data), Int32(sampling_rate), Float64(start_freq), Float64(stop_freq), Int32(order), Int32(filter_type), Float64(ripple))
    return
end

@brainflow_rethrow function remove_environmental_noise(data, sampling_rate::Integer, noise_type::EnvNoiseType)
    ccall((:remove_environmental_noise, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Cint),
            data, length(data), Int32(sampling_rate), Int32(noise_type))
    return
end

@brainflow_rethrow function perform_rolling_filter(data, period::Integer, operation::AggType)
    ccall((:perform_rolling_filter, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Cint),
            data, length(data), Int32(period), Int32(operation))
    return
end

@brainflow_rethrow function detrend(data, operation::DetType)
    ccall((:detrend, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint),
            data, length(data), Int32(operation))
    return
end

@brainflow_rethrow function perform_wavelet_denoising(data, wavelet::WaveletType, decomposition_level::Integer,
                                                      wavelet_denoising::WaveletDenoisingType,
                                                      threshold::ThresholdType,
                                                      extension::WaveletExtensionType,
                                                      noise_level::NoiseEstimationLevelType)
    ccall((:perform_wavelet_denoising, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Cint, Cint, Cint, Cint, Cint),
            data, length(data), Int32(wavelet), Int32(decomposition_level), Int32(wavelet_denoising),
            Int32(threshold), Int32(extension), Int32(noise_level))
    return
end

@brainflow_rethrow function perform_downsampling(data, period::Integer, operation::AggType)
    len = Integer(floor(length(data) / period))
    downsampled_data = Vector{Float64}(undef, len)
    ccall((:perform_downsampling, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Cint, Ptr{Float64}),
            data, length(data), Int32(period), Int32(operation), downsampled_data)
    return downsampled_data
end

@brainflow_rethrow function write_file(data, file_name::String, file_mode::String)
    shape = size(data)
    flatten = transpose(vcat(data))
    flaten = reshape(flatten, (1, shape[1] * shape[2]))
    flatten = copy(flatten)
    ccall((:write_file, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Ptr{UInt8}, Ptr{UInt8}),
            flatten, shape[1], shape[2], file_name, file_mode)
    return
end

@brainflow_rethrow function get_nearest_power_of_two(value::Integer)
    power_of_two = Vector{Cint}(undef, 1)
    ccall((:get_nearest_power_of_two, DATA_HANDLER_INTERFACE), Cint, (Cint, Ptr{Cint}),
                Int32(value), power_of_two)
    return power_of_two[1]
end

@brainflow_rethrow function calc_stddev(data)
    output = Vector{Float64}(undef, 1)
    ccall((:calc_stddev, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Ptr{Float64}),
                data, 0, length(data), output)
    return output[1]
end

@brainflow_rethrow function get_num_elements_in_file(file_name::String)
    num_elements = Vector{Cint}(undef, 1)
    ccall((:get_num_elements_in_file, DATA_HANDLER_INTERFACE), Cint, (Ptr{UInt8}, Ptr{Cint}),
                file_name, num_elements)
    return num_elements
end

@brainflow_rethrow function read_file(file_name::String)
    num_elements = get_num_elements_in_file(file_name)
    data = Vector{Float64}(undef, num_elements[1])
    num_cols = Vector{Cint}(undef, 1)
    num_rows = Vector{Cint}(undef, 1)
    ccall((:read_file, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Ptr{Cint}, Ptr{Cint}, Ptr{UInt8}, Cint),
                data, num_rows, num_cols, file_name, num_elements[1])
    value = transpose(reshape(data, (num_cols[1], num_rows[1])))
    return value
end

@brainflow_rethrow function perform_wavelet_transform(data, wavelet::WaveletType, decomposition_level::Integer, extension::WaveletExtensionType)
    wavelet_coeffs = Vector{Float64}(undef, length(data) + 2 * (40 + 1))
    lengths = Vector{Cint}(undef, decomposition_level + 1)
    ccall((:perform_wavelet_transform, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Cint, Cint, Ptr{Float64}, Ptr{Cint}),
            data, length(data), Int32(wavelet), Int32(decomposition_level), Int32(extension), wavelet_coeffs, lengths)
    return wavelet_coeffs[1:sum(lengths)], lengths
end

@brainflow_rethrow function perform_inverse_wavelet_transform(wavelet_output, original_data_len::Integer, wavelet::WaveletType, decomposition_level::Integer, extension::WaveletExtensionType)
    original_data = Vector{Float64}(undef, original_data_len)
    ccall((:perform_inverse_wavelet_transform, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Cint, Cint, Ptr{Float64}, Ptr{Float64}),
            wavelet_output[1], Int32(original_data_len), Int32(wavelet), Int32(decomposition_level), Int32(extension), wavelet_output[2], original_data)
    return original_data
end

@brainflow_rethrow function get_csp(data, labels)
    n_epochs = size(data, 1)
    n_channels = size(data, 2)
    n_times = size(data, 3)

    temp_data1d = Vector{Float64}(undef, Integer(n_epochs * n_channels * n_times))
    for e=1:n_epochs
        for c=1:n_channels
            for t=1:n_times
                temp_data1d[(e-1) * n_channels * n_times + (c-1) * n_times + t] = data[e, c, t]
            end
        end
    end

    temp_filters = Vector{Float64}(undef, Integer(n_channels * n_channels))
    output_eigenvalues = Vector{Float64}(undef, Integer(n_channels))

    ccall((:get_csp, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Ptr{Float64}, Cint, Cint, Cint, Ptr{Float64}, Ptr{Float64}), temp_data1d, labels, Int32(n_epochs), Int32(n_channels), Int32(n_times), temp_filters, output_eigenvalues)

    output_filters = Array{Float64,2}(undef, n_channels, n_channels)
    for i=1:n_channels
        for j=1:n_channels
            output_filters[i, j] = temp_filters[(i-1) * n_channels + j]
        end
    end
    return output_filters, output_eigenvalues
end

@brainflow_rethrow function get_window(window_function::WinType, window_len::Integer)
    window_data = Vector{Float64}(undef, Integer(window_len))
    ccall((:get_window, DATA_HANDLER_INTERFACE), Cint, (Cint, Cint, Ptr{Float64}),
    Int32(window_function), Int32(window_len), window_data)
    return window_data
end

@brainflow_rethrow function perform_fft(data, window::WinType)

    function is_power_of_two(value)
        (value != 0) && (value & (value - 1) == 0)
    end

    if !is_power_of_two(length(data))
        throw(BrainFlowError(string("Data Len must be power of two ", INVALID_ARGUMENTS_ERROR), Integer(INVALID_ARGUMENTS_ERROR)))
    end

    temp_re = Vector{Float64}(undef, Integer(length(data) / 2) + 1)
    temp_im = Vector{Float64}(undef, Integer(length(data) / 2) + 1)
    res = Vector{Complex}(undef, Integer(length(data) / 2) + 1)

    ccall((:perform_fft, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Ptr{Float64}, Ptr{Float64}),
            data, length(data), Int32(window), temp_re, temp_im)
    for i in 1:Integer(length(data) / 2) + 1
        res[i] = Complex(temp_re[i], temp_im[i])
    end
    return res
end

@brainflow_rethrow function perform_ifft(data)

    temp_re = Vector{Float64}(undef, length(data))
    temp_im = Vector{Float64}(undef, length(data))
    res = Vector{Float64}(undef, 2 * (length(data) - 1))

    for i in 1:length(data)
        temp_re[i] = data[i].re
        temp_im[i] = data[i].im
    end

    ccall((:perform_ifft, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Ptr{Float64}, Cint, Ptr{Float64}),
            temp_re, temp_im, length(res), res)
    return res
end

function get_avg_band_powers(data, channels, sampling_rate::Integer, apply_filter::Bool)
    bands = [(2.0, 4.0), (4.0, 8.0), (8.0, 13.0), (13.0, 30.0), (30.0, 45.0)]
    return get_custom_band_powers(data, bands, channels, sampling_rate, apply_filter)
end

@brainflow_rethrow function get_custom_band_powers(data, bands, channels, sampling_rate::Integer, apply_filter::Bool)

    shape = size(data)
    data_1d = reshape(transpose(data[channels,:]), (1, size(channels)[1] * shape[2]))
    data_1d = copy(data_1d)

    start_freqs = [first(p) for p in bands]
    stop_freqs = [last(p) for p in bands]

    temp_avgs = Vector{Float64}(undef, length(start_freqs))
    temp_stddevs = Vector{Float64}(undef, length(start_freqs))

    ccall((:get_custom_band_powers, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Ptr{Float64}, Ptr{Float64}, Cint, Cint, Cint, Ptr{Float64}, Ptr{Float64}),
            data_1d, size(channels)[1], shape[2], start_freqs, stop_freqs, length(start_freqs), Int32(sampling_rate), Int32(apply_filter), temp_avgs, temp_stddevs)
    return temp_avgs, temp_stddevs
end

@brainflow_rethrow function get_psd(data, sampling_rate::Integer, window::WinType)

    function is_power_of_two(value)
        (value != 0) && (value & (value - 1) == 0)
    end

    if !is_power_of_two(length(data))
        throw(BrainFlowError(string("Data Len must be power of two ", INVALID_ARGUMENTS_ERROR), Integer(INVALID_ARGUMENTS_ERROR)))
    end

    temp_ampls = Vector{Float64}(undef, Integer(length(data) / 2) + 1)
    temp_freqs = Vector{Float64}(undef, Integer(length(data) / 2) + 1)

    ccall((:get_psd, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Cint, Ptr{Float64}, Ptr{Float64}),
            data, length(data), Int32(sampling_rate), Int32(window), temp_ampls, temp_freqs)
    return temp_ampls, temp_freqs
end

@brainflow_rethrow function get_psd_welch(data, nfft::Integer, overlap::Integer, sampling_rate::Integer, window::WinType)

    function is_power_of_two(value)
        (value != 0) && (value & (value - 1) == 0)
    end

    if !is_power_of_two(nfft)
        throw(BrainFlowError(string("nfft must be power of two ", INVALID_ARGUMENTS_ERROR), Integer(INVALID_ARGUMENTS_ERROR)))
    end

    temp_ampls = Vector{Float64}(undef, Integer(nfft / 2) + 1)
    temp_freqs = Vector{Float64}(undef, Integer(nfft / 2) + 1)

    ccall((:get_psd_welch, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Cint, Cint, Cint, Ptr{Float64}, Ptr{Float64}),
            data, length(data), Int32(nfft), Int32(overlap), Int32(sampling_rate), Int32(window), temp_ampls, temp_freqs)
    return temp_ampls, temp_freqs
end

@brainflow_rethrow function get_band_power(psd, freq_start::Float64, freq_end::Float64)
    band_power = Vector{Float64}(undef, 1)
    ccall((:get_band_power, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Ptr{Float64}, Cint, Float64, Float64, Ptr{Float64}),
            psd[1], psd[2], length(psd[1]), Float64(freq_start), Float64(freq_end), band_power)
    return band_power[1]
end


text:
import .Luxor: Drawing, finish, Point, background, sethue, fontsize, rect, text

prof_to_svg(::Any, ::Any, ::Any; kwargs...) = nothing

combine_gantt_images(::Any, ::Any, ::Any, ::Any) = ("", "")

function draw_gantt(ctx, svg_path, prof_path; delay=2, width=1000, height=640, window_length=20)
    root_time = time_ns()
    window_logs = []
    if window_length !== nothing
        window_start = -window_length
    else
        window_start = 0
    end
    window_finish = 0
    procs = Set()
    for gp in ctx.procs
        for proc in Dagger.get_processors(gp)
            push!(procs, proc)
        end
    end
    if (height/length(procs)) < 50
        height = length(procs) * 50
        @warn "SVG height too small; resizing to $height pixels"
    end

    image_idx = 0
    while continue_rendering[]
        sleep(delay)
        logs = Dagger.get_logs!(ctx.log_sink, true) # get raw events
        update_window_logs!(window_logs, logs; root_time=root_time, window_start=window_start)
        isempty(window_logs) && continue

        # Concatenate and render profile data
        if ctx.profile
            prof_data, prof_lidict = logs_to_stackframes(window_logs)
            prof_to_svg(prof_path, prof_data, prof_lidict, image_idx; width=width)
        end

        for proc in unique(map(x->x[1].timeline[2], filter(x->x[1].category==:compute, window_logs)))
            push!(procs, proc)
        end
        colors = Colors.distinguishable_colors(length(procs))
        procs_len = length(procs)
        proc_height = height/(3procs_len)
        @debug "Gantt: Start"
        if isfile(svg_path)
            Drawing(width, height, svg_path)
        else
            Drawing(width, height, joinpath(svg_path, repr(image_idx) * ".svg"))
        end
        background("white")
        for (proc_idx, proc) in enumerate(sort(collect(procs); lt=proclt))
            ypos = (proc_idx-0.5)*(height/procs_len)
            sethue("grey15")
            fontsize(round(Int,proc_height))
            text(getname(proc), Point(width/2,ypos-(proc_height/2)))
            rect(Point(1,ypos+(proc_height/3)),width-2,proc_height,:stroke)
            fontsize(8)
            text("$(window_start) s", Point(0,ypos); halign=:left)
            text("$(window_finish) s", Point(width-8,ypos); halign=:right)
            proc_color = colors[proc_idx]
            for log in filter(x->x[1].timeline[2]==proc, filter(x->x[1].category==:compute, window_logs))
                length(log) == 1 && log[1] isa Dagger.Event{:finish} && error("Unpaired finish!")
                log_start_s = (log[1].timestamp-root_time)/(1000^3)
                log_finish_s = if length(log) == 2
                    log_finish_s = (log[2].timestamp-root_time)/(1000^3)
                else
                    window_finish+1
                end
                xstart = ((log_start_s-window_start)/(window_finish-window_start))*width
                xfinish = ((log_finish_s-window_start)/(window_finish-window_start))*width
                sethue(proc_color)
                rect(Point(xstart,ypos+(proc_height/3)+1),xfinish-xstart,proc_height-2,:fill)
                sethue("black")
                rect(Point(xstart,ypos+(proc_height/3)+1),xfinish-xstart,proc_height-2,:stroke)
            end
            for log in filter(x->x[1].category==:scheduler_init, window_logs)
                log_start_s = (log[1].timestamp-root_time)/(1000^3)
                log_finish_s = if length(log) == 2
                    log_finish_s = (log[2].timestamp-root_time)/(1000^3)
                else
                    window_finish+1
                end
                xstart = ((log_start_s-window_start)/(window_finish-window_start))*width
                xfinish = ((log_finish_s-window_start)/(window_finish-window_start))*width
                sethue("red")
                rect(Point(xstart,0),#=xfinish-xstart=#1,height)
            end
        end
        finish()
        @debug "Gantt: Finish"
        window_finish = (time_ns()-root_time)/(1000^3)
        if window_length !== nothing
            window_start = window_finish-window_length
        end
        image_idx += 1
    end
    if isdir(svg_path)
        final_paths = try
            combine_gantt_images(ctx, svg_path, prof_path, delay)
        catch err
            @error "Image-to-video failed" exception=err
            ("", "")
        end
        continue_rendering[] = true
        put!(render_results, final_paths)
    end
end

function show_gantt(ctx; delay=2, port=8000, width=1000, height=640, window_length=20, live=false)
    svg_prefix = "/tmp/dagger-gantt-$(String(rand(UInt8(97):UInt8(122), 16)))"
    prof_prefix = "/tmp/dagger-gantt-$(String(rand(UInt8(97):UInt8(122), 16)))"
    if live
        # Create SVG file
        # TODO: Ensure this gets cleaned up
        svg_path = svg_prefix*".svg"
        prof_path = prof_prefix*".svg"
        open(svg_path, "w") do io
            write(io, "
Waiting for Gantt data...
")
        end
        open(prof_path, "w") do io
            write(io, "
Waiting for profile data...
")
        end
    else
        svg_path = svg_prefix
        prof_path = prof_prefix
        mkpath(svg_path)
        mkpath(prof_path)
    end

    # Start drawing task
    Threads.@spawn begin
        try
            draw_gantt(ctx, svg_path, prof_path; delay=delay, width=width,
                       height=height, window_length=window_length)
        catch err
            @error exception=(err,catch_backtrace())
        end
    end

    if live
        serve_gantt(svg_path, prof_path; port=port, delay=delay)
    end
end

ext: jl


text:
export AuditLog

const AUDIT_LOG_CHANGE_TYPES = Dict(
    "name"                          => (String, Guild),
    "icon_hash"                     => (String, Guild) ,
    "splash_hash"                   => (String, Guild),
    "owner_id"                      => (Snowflake, Guild),
    "region"                        => (String, Guild),
    "afk_channel_id"                => (Snowflake, Guild),
    "afk_timeout"                   => (Int, Guild),
    "mfa_level"                     => (MFALevel, Guild),
    "verification_level"            => (VerificationLevel, Guild),
    "explicit_content_filter"       => (ExplicitContentFilterLevel, Guild),
    "default_message_notifications" => (MessageNotificationLevel, Guild),
    "vanity_url_code"               => (String, Guild),
    "\$add"                         => (Vector{Role}, Guild),
    "\$remove"                      => (Vector{Role}, Guild),
    "prune_delete_days"             => (Int, Guild),
    "widget_enabled"                => (Bool, Guild),
    "widget_channel_id"             => (Snowflake, Guild),
    "position"                      => (Int, DiscordChannel),
    "topic"                         => (String, DiscordChannel),
    "bitrate"                       => (Int, DiscordChannel),
    "permission_overwrites"         => (Vector{Overwrite}, DiscordChannel),
    "nsfw"                          => (Bool, DiscordChannel),
    "application_id"                => (Snowflake, DiscordChannel),
    "permissions"                   => (Int, Role),
    "color"                         => (Int, Role),
    "hoist"                         => (Bool, Role),
    "mentionable"                   => (Bool, Role),
    "allow"                         => (Int, Role),
    "deny"                          => (Int, Role),
    "code"                          => (String, Invite),
    "channel_id"                    => (Snowflake, Invite),
    "inviter_id"                    => (Snowflake, Invite),
    "max_uses"                      => (Int, Invite),
    "uses"                          => (Int, Invite),
    "max_age"                       => (Int, Invite),
    "temporary"                     => (Bool, Invite),
    "deaf"                          => (Bool, User),
    "mute"                          => (Bool, User),
    "nick"                          => (String, User),
    "avatar_hash"                   => (String, User),
    "id"                            => (Snowflake, Any),
    "type"                          => (Any, Any),
    # Undocumented.
    "rate_limit_per_user"           => (Int, DiscordChannel),
)

"""
[`AuditLog`](@ref) action types.
More details [here](https://discordapp.com/developers/docs/resources/audit-log#audit-log-entry-object-audit-log-events).
"""
@enum ActionType  begin
    AT_GUILD_UPDATE=1
    AT_CHANNEL_CREATE=10
    AT_CHANNEL_UPDATE=11
    AT_CHANNEL_DELETE=12
    AT_CHANNEL_OVERWRITE_CREATE=13
    AT_CHANNEL_OVERWRITE_UPDATE=14
    AT_CHANNEL_OVERWRITE_DELETE=15
    AT_MEMBER_KICK=20
    AT_MEMBER_PRUNE=21
    AT_MEMBER_BAN_ADD=22
    AT_MEMBER_BAN_REMOVE=23
    AT_MEMBER_UPDATE=24
    AT_MEMBER_ROLE_UPDATE=25
    AT_ROLE_CREATE=30
    AT_ROLE_UPDATE=31
    AT_ROLE_DELETE=32
    AT_INVITE_CREATE=40
    AT_INVITE_UPDATE=41
    AT_INVITE_DELETE=42
    AT_WEBHOOK_CREATE=50
    AT_WEBHOOK_UPDATE=51
    AT_WEBHOOK_DELETE=52
    AT_EMOJI_CREATE=60
    AT_EMOJI_UPDATE=61
    AT_EMOJI_DELETE=62
    AT_MESSAGE_DELETE=72
end
@boilerplate ActionType :lower

"""
A change item in an [`AuditLogEntry`](@ref).

The first type parameter is the type of `new_value` and `old_value`. The second is the type
of the entity that `new_value` and `old_value` belong(ed) to.

More details [here](https://discordapp.com/developers/docs/resources/audit-log#audit-log-change-object).
"""
struct AuditLogChange{T, U}
    new_value::Union{T, Missing}
    old_value::Union{T, Missing}
    key::String
    type::Type{U}
end
@boilerplate AuditLogChange :docs

function AuditLogChange(d::Dict{String, Any})
    return if haskey(AUDIT_LOG_CHANGE_TYPES, d["key"])
        T, U = AUDIT_LOG_CHANGE_TYPES[d["key"]]
        func = if T === Any
            identity
        elseif T === Snowflake
            snowflake
        elseif T <: Vector
            eltype(T)
        else
            T
        end

        new_value = if haskey(d, "new_value")
            d["new_value"] isa Vector ? func.(d["new_value"]) : func(d["new_value"])
        else
            missing
        end
        old_value = if haskey(d, "old_value")
            d["old_value"] isa Vector ? func.(d["old_value"]) : func(d["old_value"])
        else
            missing
        end

        AuditLogChange{T, U}(new_value, old_value, d["key"], U)
    else
        AuditLogChange{Any, Any}(
            get(d, "new_value", missing),
            get(d, "old_value", missing),
            d["key"],
            Any,
        )
    end
end

function JSON.lower(x::AuditLogChange)
    d = Dict()
    if !ismissing(x.new_value)
        d["new_value"] = x.new_value
    end
    if !ismissing(x.old_value)
        d["old_value"] = x.new_value
    end
    d["key"] = x.key
    return d
end

"""
Optional information in an [`AuditLogEntry`](@ref).
More details [here](https://discordapp.com/developers/docs/resources/audit-log#audit-log-entry-object-optional-audit-entry-info).
"""
struct AuditLogOptions
    delete_member_days::Union{Int, Missing}
    members_removed::Union{Int, Missing}
    channel_id::Union{Snowflake, Missing}
    count::Union{Int, Missing}
    id::Union{Snowflake, Missing}
    type::Union{OverwriteType, Missing}
    role_name::Union{String, Missing}
end
@boilerplate AuditLogOptions :docs :merge

function AuditLogOptions(d::Dict{String, Any})
    return AuditLogOptions(
        haskey(d, "delete_member_days") ? parse(Int, d["delete_member_days"]) : missing,
        haskey(d, "members_removed") ? parse(Int, d["members_removed"]) : missing,
        haskey(d, "channel_id") ? snowflake(d["channel_id"]) : missing,
        haskey(d, "count") ? parse(Int, d["count"]) : missing,
        haskey(d, "id") ? snowflake(d["id"]) : missing,
        haskey(d, "type") ? OverwriteType(d["type"]) : missing,
        get(d, "role_name", missing),
    )
end

function JSON.lower(x::AuditLogOptions)
    d = Dict()
    if !ismissing(x.delete_member_days)
        d["delete_member_days"] = string(d.delete_member_days)
    end
    if !ismissing(x.members_removed)
        d["members_removed"] = string(d.members_removed)
    end
    if !ismissing(x.channel_id)
        d["channel_id"] = string(d.members_removed)
    end
    if !ismissing(x.count)
        d["count"] = string(x.count)
    end
    if !ismissing(x.id)
        d["id"] = string(x.id)
    end
    if !ismissing(x.type)
        d["type"] = string(x.type)
    end
    if !ismissing(x.role_name)
        d["role_name"] = x.role_name
    end
    return d
end

"""
An entry in an [`AuditLog`](@ref).
More details [here](https://discordapp.com/developers/docs/resources/audit-log#audit-log-entry-object).
"""
struct AuditLogEntry
    target_id::Union{Snowflake, Nothing}
    changes::Union{Vector{AuditLogChange}, Missing}
    user_id::Snowflake
    id::Snowflake
    action_type::ActionType
    options::Union{AuditLogOptions, Missing}
    reason::Union{String, Missing}
end
@boilerplate AuditLogEntry :dict :docs :lower :merge

"""
An audit log.
More details [here](https://discordapp.com/developers/docs/resources/audit-log#audit-log-object).
"""
struct AuditLog
    webhooks::Vector{Webhook}
    users::Vector{User}
    audit_log_entries::Vector{AuditLogEntry}
end
@boilerplate AuditLog :dict :docs :lower :merge


text:
"""
    PrescribedConditions{T}

Describes the forces, moments, displacements, and/or rotations prescribed at a
point at a specific point in time.
"""
struct PrescribedConditions{T}
    force::SVector{6, Bool}
    value::SVector{6, T}
    follower::SVector{6, T}
end
Base.eltype(::PrescribedConditions{T}) where T = T

function PrescribedConditions{T}(p::PrescribedConditions) where T
    PrescribedConditions{T}(p.force, p.value, p.follower)
end

PrescribedConditions(force, value, follower) = PrescribedConditions(force, promote(value, follower)...)

"""
    PrescribedConditions(t = 0.0; kwargs...)

Return the combined prescribed conditions at a point (as an object of type
[`PrescribedConditions`](@ref)) at time `t`.

Individual prescribed conditions may be assigned as either a scalar parameter or
as a function of time.

Prescribed Wiener-Milenkovic parameters must satisfy the following inequality:
sqrt(theta_x^2 + theta_y^2 + theta_z^2) <= 4.  Note that this restriction still
allows all possible rotations to be represented.

# Keyword Arguments
 - `ux`: Prescribed x-direction displacement of the point
 - `uy`: Prescribed y-direction displacement of the point
 - `uz`: Prescribed z-direction displacement of the point
 - `theta_x`: Prescribed first Wiener-Milenkovic parameter of the point
 - `theta_y`: Prescribed second Wiener-Milenkovic parameter of the point
 - `theta_z`: Prescribed third Wiener-Milenkovic parameter of the point
 - `Fx`: Prescribed force in x-direction applied on the point
 - `Fy`: Prescribed force in y-direction applied on the point
 - `Fz`: Prescribed force in z-direction applied on the point
 - `Mx`: Prescribed moment about x-axis applied on the point
 - `My`: Prescribed moment about y-axis applied on the point
 - `Mz`: Prescribed moment about z-axis applied on the point
 - `Fx_follower`: Prescribed follower force in x-direction applied on the point
 - `Fy_follower`: Prescribed follower force in y-direction applied on the point
 - `Fz_follower`: Prescribed follower force in z-direction applied on the point
 - `Mx_follower`: Prescribed follower moment about x-axis applied on the point
 - `My_follower`: Prescribed follower moment about y-axis applied on the point
 - `Mz_follower`: Prescribed follower moment about z-axis applied on the point
"""
function PrescribedConditions(t = 0.0;
    ux = nothing,
    uy = nothing,
    uz = nothing,
    theta_x = nothing,
    theta_y = nothing,
    theta_z = nothing,
    Fx = nothing,
    Fy = nothing,
    Fz = nothing,
    Mx = nothing,
    My = nothing,
    Mz = nothing,
    Fx_follower = nothing,
    Fy_follower = nothing,
    Fz_follower = nothing,
    Mx_follower = nothing,
    My_follower = nothing,
    Mz_follower = nothing)

    # Set function for first slot
    if isnothing(ux)
        f1 = true
        v1 = isnothing(Fx) ? 0.0 : typeof(Fx) <: Number ? Fx : Fx(t)
        v1_f = isnothing(Fx_follower) ? 0.0 : typeof(Fx_follower) <: Number ? Fx_follower : Fx_follower(t)
    elseif !isnothing(ux) && (isnothing(Fx) && isnothing(Fx_follower))
        f1 = false
        v1 = isnothing(ux) ? 0.0 : typeof(ux) <: Number ? ux : ux(t)
        v1_f = 0.0
    else
        error("Both `ux` and `Fx` or `Fx_follower` cannot be specified at the same time")
    end

    # Set function for second slot
    if isnothing(uy)
        f2 = true
        v2 = isnothing(Fy) ? 0.0 : typeof(Fy) <: Number ? Fy : Fy(t)
        v2_f = isnothing(Fy_follower) ? 0.0 : typeof(Fy_follower) <: Number ? Fy_follower : Fy_follower(t)
    elseif !isnothing(uy) && (isnothing(Fy) && isnothing(Fy_follower))
        f2 = false
        v2 = isnothing(uy) ? 0.0 : typeof(uy) <: Number ? uy : uy(t)
        v2_f = 0.0
    else
        error("Both `uy` and `Fy` or `Fy_follower` cannot be specified at the same time")
    end

    # Set function for third slot
    if isnothing(uz)
        f3 = true
        v3 = isnothing(Fz) ? 0.0 : typeof(Fz) <: Number ? Fz : Fz(t)
        v3_f = isnothing(Fz_follower) ? 0.0 : typeof(Fz_follower) <: Number ? Fz_follower : Fz_follower(t)
    elseif !isnothing(uz) && (isnothing(Fz) && isnothing(Fz_follower))
        f3 = false
        v3 = isnothing(uz) ? 0.0 : typeof(uz) <: Number ? uz : uz(t)
        v3_f = 0.0
    else
        error("Both `uz` and `Fz` or `Fz_follower` cannot be specified at the same time")
    end

    # Set function for fourth slot
    if isnothing(theta_x)
        f4 = true
        v4 = isnothing(Mx) ? 0.0 : typeof(Mx) <: Number ? Mx : Mx(t)
        v4_f = isnothing(Mx_follower) ? 0.0 : typeof(Mx_follower) <: Number ? Mx_follower : Mx_follower(t)
    elseif !isnothing(theta_x) && (isnothing(Mx) && isnothing(Mx_follower))
        f4 = false
        v4 = isnothing(theta_x) ? 0.0 : typeof(theta_x) <: Number ? theta_x : theta_x(t)
        v4_f = 0.0
    else
        error("Both `theta_x` and `Mx` or `Mx_follower` cannot be specified at the same time")
    end

    # Set function for fifth slot
    if isnothing(theta_y)
        f5 = true
        v5 = isnothing(My) ? 0.0 : typeof(My) <: Number ? My : My(t)
        v5_f = isnothing(My_follower) ? 0.0 : typeof(My_follower) <: Number ? My_follower : My_follower(t)
    elseif !isnothing(theta_y) && (isnothing(My) && isnothing(My_follower))
        f5 = false
        v5 = isnothing(theta_y) ? 0.0 : typeof(theta_y) <: Number ? theta_y : theta_y(t)
        v5_f = 0.0
    else
        error("Both `theta_y` and `My` or `My_follower` cannot be specified at the same time")
    end

    # Set function for sixth slot
    if isnothing(theta_z)
        f6 = true
        v6 = isnothing(Mz) ? 0.0 : typeof(Mz) <: Number ? Mz : Mz(t)
        v6_f = isnothing(Mz_follower) ? 0.0 : typeof(Mz_follower) <: Number ? Mz_follower : Mz_follower(t)
    elseif !isnothing(theta_z) && (isnothing(Mz) && isnothing(Mz_follower))
        f6 = false
        v6 = isnothing(theta_z) ? 0.0 : typeof(theta_z) <: Number ? theta_z : theta_z(t)
        v6_f = 0.0
    else
        error("Both `theta_z` and `Mz` or `Mz_follower` cannot be specified at the same time")
    end

    force = SVector(f1, f2, f3, f4, f5, f6)
    value = SVector(v1, v2, v3, v4, v5, v6)
    follower = SVector(v1_f, v2_f, v3_f, v4_f, v5_f, v6_f)

    return PrescribedConditions(force, value, follower)
end

"""
    DistributedLoads{T}

Contains the integrated distributed forces and moments for each beam element at
a specific time step.

# Fields
 - f1: Integrated non-follower distributed force for the beam element's left endpoint for each time step
 - f2: Integrated non-follower distributed force for the beam element's right endpoint for each time step
 - m1: Integrated non-follower distributed moment for the beam element's left endpoint for each time step
 - m2: Integrated non-follower distributed moment for the beam element's right endpoint for each time step
 - f1_follower: Integrated follower distributed force for the beam element's left endpoint for each time step
 - f2_follower: Integrated follower distributed force for the beam element's right endpoint for each time step
 - m1_follower: Integrated follower distributed moment for the beam element's left endpoint for each time step
 - m2_follower: Integrated follower distributed moment for the beam element's right endpoint for each time step
"""
struct DistributedLoads{T}
    f1::SVector{3, T}
    f2::SVector{3, T}
    m1::SVector{3, T}
    m2::SVector{3, T}
    f1_follower::SVector{3, T}
    f2_follower::SVector{3, T}
    m1_follower::SVector{3, T}
    m2_follower::SVector{3, T}
end
Base.eltype(::DistributedLoads{T}) where T = T

function DistributedLoads{T}(d::DistributedLoads) where T
    return DistributedLoads{T}(d.f1, d.f2, d.m1, d.m2, d.f1_follower, d.f2_follower,
        d.m1_follower, d.m2_follower)
end

DistributedLoads(f1, f2, m1, m2, f1_follower, f2_follower, m1_follower, m2_follower) =
    DistributedLoads(promote(f1, f2, m1, m2, f1_follower, f2_follower, m1_follower, m2_follower)...)

"""
    DistributedLoads(assembly, ielem; kwargs...)

Return the integrated distributed loads at a point (as an object of type [`DistributedLoads`](@ref)).

# Arguments
 - `assembly`: The beam element assembly
 - `ielem`: The index of the beam element which the distributed load is assigned to

# Keyword Arguments
 - `s1 = 0.0`: Start of the beam element (used solely for integrating the distributed loads)
 - `s2 = 1.0`: End of the beam element (used solely for integrating the distributed loads)
 - `method = (f, a, b) -> gauss_quadrature(f, a, b)`: Method which integrates function
    `f` from `a` to `b`. Defaults to the Gauss-Legendre quadrature with 4 points on each element.
 - `fx = (s) -> 0.0`: Distributed non-follower force on beam element in x-direction
 - `fy = (s) -> 0.0`: Distributed non-follower force on beam element in y-direction
 - `fz = (s) -> 0.0`: Distributed non-follower force on beam element in z-direction
 - `mx = (s) -> 0.0`: Distributed non-follower moment on beam element in x-direction
 - `my = (s) -> 0.0`: Distributed non-follower moment on beam element in y-direction
 - `mz = (s) -> 0.0`: Distributed non-follower moment on beam element in z-direction
 - `fx_follower = (s) -> 0.0`: Distributed follower force on beam element in x-direction
 - `fy_follower = (s) -> 0.0`: Distributed follower force on beam element in y-direction
 - `fz_follower = (s) -> 0.0`: Distributed follower force on beam element in z-direction
 - `mx_follower = (s) -> 0.0`: Distributed follower moment on beam element in x-direction
 - `my_follower = (s) -> 0.0`: Distributed follower moment on beam element in y-direction
 - `mz_follower = (s) -> 0.0`: Distributed follower moment on beam element in z-direction
"""
function DistributedLoads(assembly, ielem;
    s1 = 0.0,
    s2 = 1.0,
    method = gauss_quadrature,
    fx = (s) -> 0.0,
    fy = (s) -> 0.0,
    fz = (s) -> 0.0,
    mx = (s) -> 0.0,
    my = (s) -> 0.0,
    mz = (s) -> 0.0,
    fx_follower = (s) -> 0.0,
    fy_follower = (s) -> 0.0,
    fz_follower = (s) -> 0.0,
    mx_follower = (s) -> 0.0,
    my_follower = (s) -> 0.0,
    mz_follower = (s) -> 0.0,
    )

    return DistributedLoads(assembly, ielem, 0.0;
        s1 = s1,
        s2 = s2,
        method = method,
        fx = (s, t) -> fx(s),
        fy = (s, t) -> fy(s),
        fz = (s, t) -> fz(s),
        mx = (s, t) -> mx(s),
        my = (s, t) -> my(s),
        mz = (s, t) -> mz(s),
        fx_follower = (s, t) -> fx_follower(s),
        fy_follower = (s, t) -> fy_follower(s),
        fz_follower = (s, t) -> fz_follower(s),
        mx_follower = (s, t) -> mx_follower(s),
        my_follower = (s, t) -> my_follower(s),
        mz_follower = (s, t) -> mz_follower(s),
        )
end

"""
    DistrubutedLoads(assembly, ielem, t; kwargs...)

Return the integrated distributed loads on a beam element (as an object of type
[`DistributedLoads`](@ref)) at time `t`.

# Arguments
 - `assembly`: The beam element assembly
 - `ielem`: The index of the beam element which the distributed load is assigned to
 - `t`: time at which to evaluate the distributed loads

# Keyword Arguments
 - `s1 = 0.0`: Start of the beam element (used solely for integrating the distributed loads)
 - `s2 = 1.0`: End of the beam element (used solely for integrating the distributed loads)
 - `method = (f, a, b) -> gauss_quadrature(f, a, b)`: Method which integrates function
    `f` from `a` to `b`. Defaults to the Gauss-Legendre quadrature with 4 points on each element.
 - `fx = (s, t) -> 0.0`: Distributed non-follower force on beam element in x-direction
 - `fy = (s, t) -> 0.0`: Distributed non-follower force on beam element in y-direction
 - `fz = (s, t) -> 0.0`: Distributed non-follower force on beam element in z-direction
 - `mx = (s, t) -> 0.0`: Distributed non-follower moment on beam element in x-direction
 - `my = (s, t) -> 0.0`: Distributed non-follower moment on beam element in y-direction
 - `mz = (s, t) -> 0.0`: Distributed non-follower moment on beam element in z-direction
 - `fx_follower = (s, t) -> 0.0`: Distributed follower force on beam element in x-direction
 - `fy_follower = (s, t) -> 0.0`: Distributed follower force on beam element in y-direction
 - `fz_follower = (s, t) -> 0.0`: Distributed follower force on beam element in z-direction
 - `mx_follower = (s, t) -> 0.0`: Distributed follower moment on beam element in x-direction
 - `my_follower = (s, t) -> 0.0`: Distributed follower moment on beam element in y-direction
 - `mz_follower = (s, t) -> 0.0`: Distributed follower moment on beam element in z-direction
"""
function DistributedLoads(assembly, ielem, t;
    s1 = 0.0,
    s2 = 1.0,
    method = gauss_quadrature,
    fx = (s, t) -> 0.0,
    fy = (s, t) -> 0.0,
    fz = (s, t) -> 0.0,
    mx = (s, t) -> 0.0,
    my = (s, t) -> 0.0,
    mz = (s, t) -> 0.0,
    fx_follower = (s, t) -> 0.0,
    fy_follower = (s, t) -> 0.0,
    fz_follower = (s, t) -> 0.0,
    mx_follower = (s, t) -> 0.0,
    my_follower = (s, t) -> 0.0,
    mz_follower = (s, t) -> 0.0,
    )

    # get element length
    ΔL = assembly.elements[ielem].L

    # create function for general coordinate
    ξ = (s) -> (s-s1)/(s2-s1)

    # element length scaling
    scaling = ΔL/(s2-s1)

    # integrate to get f1 at this point in time
    f1 = SVector(
        method((s)->(1-ξ(s))*fx(s, t)*scaling, s1, s2),
        method((s)->(1-ξ(s))*fy(s, t)*scaling, s1, s2),
        method((s)->(1-ξ(s))*fz(s, t)*scaling, s1, s2)
    )

    # integrate to get f2 at this point in time
    f2 = SVector(
        method((s)->ξ(s)*fx(s, t)*scaling, s1, s2),
        method((s)->ξ(s)*fy(s, t)*scaling, s1, s2),
        method((s)->ξ(s)*fz(s, t)*scaling, s1, s2)
    )

    # integrate to get m1 at this point in time
    m1 = SVector(
        method((s)->(1-ξ(s))*mx(s, t)*scaling, s1, s2),
        method((s)->(1-ξ(s))*my(s, t)*scaling, s1, s2),
        method((s)->(1-ξ(s))*mz(s, t)*scaling, s1, s2)
    )

    # integrate to get m2 at this point in time
    m2 = SVector(
        method((s)->ξ(s)*mx(s, t)*scaling, s1, s2),
        method((s)->ξ(s)*my(s, t)*scaling, s1, s2),
        method((s)->ξ(s)*mz(s, t)*scaling, s1, s2)
    )

    # integrate to get f1_follower at this point in time
    f1_follower =SVector(
        method((s)->(1-ξ(s))*fx_follower(s, t)*scaling, s1, s2),
        method((s)->(1-ξ(s))*fy_follower(s, t)*scaling, s1, s2),
        method((s)->(1-ξ(s))*fz_follower(s, t)*scaling, s1, s2)
    )

    # integrate to get f2_follower at this point in time
    f2_follower = SVector(
        method((s)->ξ(s)*fx_follower(s, t)*scaling, s1, s2),
        method((s)->ξ(s)*fy_follower(s, t)*scaling, s1, s2),
        method((s)->ξ(s)*fz_follower(s, t)*scaling, s1, s2)
    )

    # integrate to get m1_follower at this point in time
    m1_follower = SVector(
        method((s)->(1-ξ(s))*mx_follower(s, t)*scaling, s1, s2),
        method((s)->(1-ξ(s))*my_follower(s, t)*scaling, s1, s2),
        method((s)->(1-ξ(s))*mz_follower(s, t)*scaling, s1, s2)
    )

    # integrate to get m2_follower at this point in time
    m2_follower = SVector(
        method((s)->ξ(s)*mx_follower(s, t)*scaling, s1, s2),
        method((s)->ξ(s)*my_follower(s, t)*scaling, s1, s2),
        method((s)->ξ(s)*mz_follower(s, t)*scaling, s1, s2)
    )

    return DistributedLoads(f1, f2, m1, m2, f1_follower, f2_follower, m1_follower, m2_follower)
end

"""
    combine_loads(l1, l2)

Combine the distributed loads in `l1` and `l2`
"""
function combine_loads(l1, l2)
    f1 = l1.f1 + l2.f1
    f2 = l1.f2 + l2.f2
    m1 = l1.m1 + l2.m1
    m2 = l1.m2 + l2.m2
    f1_follower = l1.f1_follower + l2.f1_follower
    f2_follower = l1.f2_follower + l2.f2_follower
    m1_follower = l1.m1_follower + l2.m1_follower
    m2_follower = l1.m2_follower + l2.m2_follower
    return DistributedLoads(f1, f2, m1, m2, f1_follower, f2_follower, m1_follower, m2_follower)
end


text:
should_precompile = true

# Don't edit the following! Instead change the script for `snoop_bot`.
ismultios = true
ismultiversion = true
# precompile_enclosure
@static if !should_precompile
    # nothing
elseif !ismultios && !ismultiversion
    @static if (isfile("../deps/SnoopCompile/precompile/precompile_MatLang.jl"))
        include("../deps/SnoopCompile/precompile/precompile_MatLang.jl")
        _precompile_()
    end
else
    @static if Sys.islinux()
        @static if v"1.0.0-DEV" <= VERSION <= v"1.0.9"
            @static if (isfile("../deps/SnoopCompile/precompile/linux/1.0/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/linux/1.0/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.1.0-DEV" <= VERSION <= v"1.1.9"
            @static if (isfile("../deps/SnoopCompile/precompile/linux/1.1/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/linux/1.1/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.2.0-DEV" <= VERSION <= v"1.2.9"
            @static if (isfile("../deps/SnoopCompile/precompile/linux/1.2/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/linux/1.2/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.3.0-DEV" <= VERSION <= v"1.3.9"
            @static if (isfile("../deps/SnoopCompile/precompile/linux/1.3/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/linux/1.3/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.4.0-DEV" <= VERSION <= v"1.4.9"
            @static if (isfile("../deps/SnoopCompile/precompile/linux/1.4/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/linux/1.4/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.5.0-DEV" <= VERSION <= v"1.5.9"
            @static if (isfile("../deps/SnoopCompile/precompile/linux/1.5/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/linux/1.5/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.6.0-DEV" <= VERSION <= v"1.6.9"
            @static if (isfile("../deps/SnoopCompile/precompile/linux/1.6/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/linux/1.6/precompile_MatLang.jl")
                _precompile_()
            end
        else
        end

    elseif Sys.iswindows()
        @static if v"1.0.0-DEV" <= VERSION <= v"1.0.9"
            @static if (isfile("../deps/SnoopCompile/precompile/windows/1.0/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/windows/1.0/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.1.0-DEV" <= VERSION <= v"1.1.9"
            @static if (isfile("../deps/SnoopCompile/precompile/windows/1.1/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/windows/1.1/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.2.0-DEV" <= VERSION <= v"1.2.9"
            @static if (isfile("../deps/SnoopCompile/precompile/windows/1.2/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/windows/1.2/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.3.0-DEV" <= VERSION <= v"1.3.9"
            @static if (isfile("../deps/SnoopCompile/precompile/windows/1.3/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/windows/1.3/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.4.0-DEV" <= VERSION <= v"1.4.9"
            @static if (isfile("../deps/SnoopCompile/precompile/windows/1.4/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/windows/1.4/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.5.0-DEV" <= VERSION <= v"1.5.9"
            @static if (isfile("../deps/SnoopCompile/precompile/windows/1.5/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/windows/1.5/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.6.0-DEV" <= VERSION <= v"1.6.9"
            @static if (isfile("../deps/SnoopCompile/precompile/windows/1.6/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/windows/1.6/precompile_MatLang.jl")
                _precompile_()
            end
        else
        end

    elseif Sys.isapple()
        @static if v"1.0.0-DEV" <= VERSION <= v"1.0.9"
            @static if (isfile("../deps/SnoopCompile/precompile/apple/1.0/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/apple/1.0/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.1.0-DEV" <= VERSION <= v"1.1.9"
            @static if (isfile("../deps/SnoopCompile/precompile/apple/1.1/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/apple/1.1/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.2.0-DEV" <= VERSION <= v"1.2.9"
            @static if (isfile("../deps/SnoopCompile/precompile/apple/1.2/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/apple/1.2/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.3.0-DEV" <= VERSION <= v"1.3.9"
            @static if (isfile("../deps/SnoopCompile/precompile/apple/1.3/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/apple/1.3/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.4.0-DEV" <= VERSION <= v"1.4.9"
            @static if (isfile("../deps/SnoopCompile/precompile/apple/1.4/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/apple/1.4/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.5.0-DEV" <= VERSION <= v"1.5.9"
            @static if (isfile("../deps/SnoopCompile/precompile/apple/1.5/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/apple/1.5/precompile_MatLang.jl")
                _precompile_()
            end
        elseif v"1.6.0-DEV" <= VERSION <= v"1.6.9"
            @static if (isfile("../deps/SnoopCompile/precompile/apple/1.6/precompile_MatLang.jl"))
                include("../deps/SnoopCompile/precompile/apple/1.6/precompile_MatLang.jl")
                _precompile_()
            end
        else
        end

    else
    end

end # precompile_enclosure

max_stars_repo_name: juliamatlab/MatLang

max_issues_repo_name: juliamatlab/MatLang

max_forks_repo_name: juliamatlab/MatLang


text:
using Random

using Bijectors
using Turing
using Turing: Variational

using StatsFuns: softplus, invsoftplus
f = softplus
f⁻¹ = invsoftplus

# setup for plotting
using Plots, StatsPlots, LaTeXStrings
pyplot()

# used to compute closed form expression of posterior
using ConjugatePriors

# define Turing model
@model model(x) = begin
    s ~ InverseGamma(2, 3)
    m ~ Normal(0.0, sqrt(s))  # `Normal(μ, σ)` has mean μ and variance σ², i.e. parametrize with std. not variance

    for i = 1:length(x)
        x[i] ~ Normal(m, sqrt(s))
    end
end

const seeds = [125, 245, 3]
const ad_modes = [:forward_diff, :reverse_diff]

for seed ∈ seeds
    @info seed

    for ad_mode ∈ ad_modes
        @info ad_mode
        setadbackend(ad_mode)

        # set random seed
        Random.seed!(seed)

        # generate data
        x = randn(1, 2000);

        # construct model
        m = model(x)

        # ADVI
        opt = Variational.TruncatedADAGrad()   # optimizer

        advi = ADVI(10, 100)                   # <: VariationalInference
        q = Variational.meanfield(m)           # => <: VariationalPosterior

        elbo = Variational.ELBO()              # <: VariationalObjective

        μ, σs = params(q)
        θ = vcat(μ, f⁻¹.(σs))

        history = [elbo(advi, q, m, 1000)]

        # construct animation
        anim = @animate for j = 1:25
            # global q
            Variational.optimize!(elbo, advi, q, m, θ; optimizer = opt)
            μ, ω = θ[1:length(q)], θ[length(q) + 1:end]

            q = Variational.update(q, μ, f.(ω))
            samples = rand(q, 2000)

            # quick check
            println([mean(samples, dims=2), [var(x), mean(x)]])

            # plotting code assumes (samples, dim) shape so we just transpose
            samples = transpose(samples)

            # closed form computation
            # notation mapping has been verified by explicitly computing expressions
            # in "Conjugate Bayesian analysis of the Gaussian distribution" by Murphy
            μ₀ = 0.0 # => μ
            κ₀ = 1.0 # => ν, which scales the precision of the Normal
            α₀ = 2.0 # => "shape"
            β₀ = 3.0 # => "rate", which is 1 / θ, where θ is "scale"

            # prior
            pri = NormalGamma(μ₀, κ₀, α₀, β₀)

            # posterior
            post = posterior(pri, Normal, x)

            # marginal distribution of τ = 1 / σ²
            # Eq. (90) in "Conjugate Bayesian analysis of the Gaussian distribution" by Murphy
            # `scale(post)` = θ
            p_τ = Gamma(post.shape, scale(post))
            p_σ²_pdf = z -> pdf(p_τ, 1 / z) # τ => 1 / σ²

            # marginal of μ
            # Eq. (91) in "Conjugate Bayesian analysis of the Gaussian distribution" by Murphy
            p_μ = TDist(2 * post.shape)

            μₙ = post.mu    # μ → μ
            κₙ = post.nu    # κ → ν
            αₙ = post.shape # α → shape
            βₙ = post.rate  # β → rate

            # numerically more stable but doesn't seem to have effect; issue is probably internal to
            # `pdf` which needs to compute ≈ Γ(1000)
            p_μ_pdf = z -> exp(logpdf(p_μ, (z - μₙ) * exp(- 0.5 * log(βₙ) + 0.5 * log(αₙ) + 0.5 * log(κₙ))))
            # p_μ_pdf1 = z -> pdf(p_μ, (z - μₙ) / √(βₙ / (αₙ * κₙ)))

            # posterior plots
            p1 = plot();
            density!(samples[:, 1], label = "s (ADVI)", color = :blue, linestyle = :dash);
            histogram!(samples[:, 1], label = "", normed = true, alpha = 0.3, color = :blue);

            # normalize using Riemann approx. because of (almost certainly) numerical issues
            Δ = 0.001
            r = 0.75:0.001:1.50
            norm_const = sum(p_σ²_pdf.(r) .* Δ)
            plot!(r, p_σ²_pdf, label = "s (posterior)", color = :red);
            vline!([var(x)], label = "s (data)", linewidth = 1.5, color = :black, alpha = 0.7);
            xlims!(0.5, 1.5);
            title!("$(j * advi.max_iters) steps");

            p2 = plot();
            density!(samples[:, 2], label = "m (ADVI)", color = :blue, linestyle = :dash);
            histogram!(samples[:, 2], label = "", normed = true, alpha = 0.3, color = :blue);

            # normalize using Riemann approx. because of (almost certainly) numerical issues
            Δ = 0.0001
            r = -0.1 + mean(x):Δ:0.1 + mean(x)
            norm_const = sum(p_μ_pdf.(r) .* Δ)
            plot!(r, z -> p_μ_pdf(z) / norm_const, label = "m (posterior)", color = :red);
            vline!([mean(x)], label = "m (data)", linewidth = 1.5, color = :black, alpha = 0.7);

            xlims!(-0.25, 0.25);

            # visualize evolution of objective wrt. optimization iterations
            obj = elbo(advi, q, m, 1000)
            @info "ELBO" obj
            push!(history, obj)
            p3 = plot();
            plot!(1:advi.max_iters:length(history) * advi.max_iters, history, label = "")
            title!("ELBO = $obj")

            # plot the latest 25 objective evaluations to visualize trend
            p4 = plot();
            plot!(history[max(1, end - 10):end], label = "")

            p = plot(p1, p2, p3, p4; layout = (4, 1))

            @info "[$j] Done"
            p
        end
        gif(anim, "advi_w_elbo_fps15_$(seed)_$(ad_mode).gif", fps = 15)
    end
end


text:
module Demo2

export f, f2
# package code goes here

f(x) = 2x

function f2()
	println("victory")
end

end # module


text:
# # Simple Wright-Fisher

# We can create and run simple Wright-Fisher simulations with EvoDynamics.jl. To that end, we define a single haploid species, in an unstructured space, with two single genes affecting biotic and abiotic traits, respectively.

using EvoDynamics

# A simple one-species model with no spatial structure. Model parameters are in a YAML file as follows:

# ```yml
# species:
#   1:
#     name: a
#     number of genes: 2
#     number of phenotypes: 2
#     abiotic phenotypes: [1]
#     biotic phenotypes: [2]
#     migration phenotype: 0
#     migration threshold: 3.5
#     vision radius: 0
#     check fraction: 0
#     ploidy: 1
#     epistasis matrix: [1.0, 0.0, 0.0, 1.0]
#     pleiotropy matrix: [1, 0, 0, 1]
#     growth rate: 1.0
#     expression array: [0.28, 0.46]
#     selection coefficient: 0.5
#     mutation probabilities: [0.9, 0.0, 0.0]
#     mutation magnitudes: [0.05, 0.0, 0.0]
#     N: [100]
#     environmental noise: 0.01
#     optimal phenotype values:
#       - [1.76]
#     optimal phenotypes: [1, 1, 1, 1, 1, 1]
#     age: 2
#     recombination: 0
#     initial energy: 0

# model:
#   generations: 5
#   space: [1,1]
#   metric: chebyshev
#   periodic: false
#   resources: [200]
#   interactions: [-0.1]
#   food sources: [1.0]
#   seed: Null
# ```

param_file = "../../examples/paramfile1.yml" #hide
agentdata, modeldata, model = runmodel(param_file);

modeldata


text:
using CardinalDicts
using Base.Test
import CardinalDicts: pair, unpair

#=
  create a CardinalDict
      with indices 1:20 that holds Int64 values
  check length
  confirm keymax
  populate it
  use it
  unset an index
  check it
  reassign an indexable value
=#

factorials = CardinalDict{Int64}(20);

@test length(factorials) == 0
@test keymax(factorials) == 20
@test isempty(factorials) == true

for i in 1:20
    setindex!(factorials, factorial(i), i)
end

@test length(factorials) == 20
@test keymax(factorials) == 20

@test haskey(factorials, 17) == true
@test factorials[17] == factorial(17)

delete!(factorials, 17)
@test haskey(factorials, 17) == false
@test get(factorials, 17, 0) == 0

factorials[17] = factorial(17)
@test factorials[17] == factorial(17)

vec = [1.0, 3.0, 2.0]
fromvec = CardinalDict(vec)
@test fromvec[2] == 3.0
@test keymax(fromvec) == 3
@test isfull(fromvec) == true

# Container, Associative interface

tenfold = CardinalDict{String}(40)
@test length(tenfold) == 0
@test endof(tenfold) == 0
@test keymax(tenfold) == 40
@test keys(tenfold) == []
@test values(tenfold) == []

tenfold[20] = "200"
tenfold[25] = "250"
tenfold[26] = "260"

@test length(tenfold) == 3
@test endof(tenfold) == 26
@test keymax(tenfold) == 40
@test keys(tenfold) == Int8[20, 25, 26]
@test values(tenfold) == String["200", "250", "260"]
@test eltype(tenfold) == Pair{Int8, String}

delete!(tenfold, 20)
@test haskey(tenfold, 20) == false
@test get(tenfold, 20, "0") == "0"

tenfold[20] = "200"
@test haskey(tenfold, 20) == true
@test get(tenfold, 20, "0") == "200"

# etc
#=
tenfold2 = eval(parse(string(tenfold)))
@test tenfold == tenfold2
=#

#=
  define type Stringy that holds a string
  create a CardinalPairDict
      with paired indices (1:5, 1:4) that holds Stringy values
  check length
  confirm keymax
  populate it
  use it
  unset an index
  check it
  reassign an indexable value
=#

mutable struct Stringy
    value::String
end
value(x::Stringy) = x.value
function value(x::Stringy, s::String)
    x.value = s
    return x
end
Base.:(==)(x::Stringy, y::Stringy) = value(x) == value(y)

nrows = 5; ncols = 4; n = ncols*nrows
ints = collect(1:n)
matrix_of_ints = reshape(ints, nrows, ncols)
matrix_of_strs = string.(matrix_of_ints)

pairdict = CardinalPairDict{Stringy}(size(matrix_of_strs)...);

@test length(pairdict) == 0
@test keymax(pairdict) == pair(nrows, ncols)
@test isempty(pairdict) == true
@test isfull(pairdict)  == false

for r in 1:nrows
    for c in 1:ncols
        stringy = Stringy(getindex(matrix_of_strs, r, c))
        setindex!(pairdict, stringy, r, c)
    end
end

@test length(pairdict) == n
@test keymax(pairdict) == pair(nrows, ncols)
@test isempty(pairdict) == false
@test isfull(pairdict)  == false

@test haskey(pairdict, 3, 2) == true
@test value(getindex(pairdict, 3, 2)) == matrix_of_strs[3,2]

frommat = CardinalPairDict(Stringy.(matrix_of_strs))
@test keys(frommat) == keys(pairdict)
@test values(frommat) == values(pairdict)

delete!(pairdict, 3, 2)
@test haskey(pairdict, 3, 2) == false
@test get(pairdict, 3, 2, Stringy("0")) == Stringy("0")


text:
using RadiativeTransfer
using RadiativeTransfer.CrossSection

##
## STEP 1: Get the Hitran data
##

# If you have a fixed-width file stored locally
hitran_data = read_hitran("path/to/file", mol=2, iso=1, ν_min=6000, ν_max=6400)

# If you would like to download the file from the hitran database
hitran_data = read_hitran(artifact"hitran_molec_id_2_CO2", mol=2, iso=1, ν_min=6000, ν_max=6400)

##
## STEP 2: Create a model from parameters
##

# These are some example models, but you can create/customize your model however you'd like.
# Please make_hitran_model documentation for optional arguments
model_doppler = make_hitran_model(hitran_data, Doppler())
model_lorentz = make_hitran_model(hitran_data, Lorentz())
model_voigt_CPU = make_hitran_model(hitran_data, Voigt())
model_voigt_GPU = make_hitran_model(hitran_data, Voigt(), architecture=Architectures.GPU())

# If you would prefer to create an interpolation, and then interpolate the cross-section
# at other wavelengths/pressures/temperatures, you can use make_interpolation_model as such:

ν_grid = 6000:0.01:6400
pressures = 250:250:1250
temperatures = 100:75:400

model_interp = make_interpolation_model(hitran_data, Voigt(), ν_grid, pressures, temperatures)

# Note: please see make_interpolation_model documentation for optional parameters

##
## STEP 3: Calculate the absorption cross section with the created model
##

# You can specify the wavelength grid, pressure, and temperature.
absorption_cross_section(model_*, 6000:0.01:6400, 1000.1, 296.1)


text:
module pure

export
    split_backslash_last,
    generate_filepath,
    prefix_filename,
    prefix_artist_name,
    find_number_map,
    print_multiple_imgs,
    url_given_size,
    post_title,
    medium_urls,
    post_titles_in_page,
    page_urls_in_post,
    change_url_to_full,
    process_user_url,
    process_artwork_url

function split_backslash_last(string)
    return split(string, '/')[end]
end

function generate_filepath(filename)
    return "$(homedir())/Downloads/$filename"
end

function prefix_filename(old_name, new_name, number)
    img_ext = split(old_name, ".")[end]
    number_prefix = lpad(number, 3, "0")
    new_file_name = "$(number_prefix)_$(new_name).$img_ext"
    return new_file_name
end

function prefix_artist_name(name, number)
    number_prefix = lpad(number, 2, "0")
    new_file_name = "$number_prefix\n$(" " ^ 19)$name"
    return new_file_name
end

function find_number_map(x, y)
    """Translates 1-based-index coordinates into (0-) indexable number
    5 columns and 6 rows == 30 images
    -1 accounts for the input being 1-based-index but python being 0-based
    mod 5: x is cyclic for every 5
    +5y: adding a 5 for every row 'moves one square down' on the 5x6 grid

    >>> a = [find_number_map(x,y) for y in range(1,7) for x in range(1,6)]
    >>> assert a == list(range(30))
    """
    # In brackets: the accepted ranges for 5x6 grid
    if !(1 <= x <= 5 && 1 <= y <= 6)
        println("Invalid number!")
        return false
    end
    return ((x - 1) % 5) + (5 * (y - 1))
end

function print_multiple_imgs(illusts_json)
    for (index, post) in enumerate(illusts_json)
        pages = post["page_count"]
        if pages > 1
            print("#$index has $pages pages, ")
        end
    end
    println("")
end

function url_given_size(post_json, size)
    """
    size : str
        One of: ("square-medium", "medium", "large")
    """
    return post_json["image_urls"][size]
end

function post_title(current_page_illusts, post_number)
    return current_page_illusts[post_number]["title"]
end

function medium_urls(current_page_illusts)
    urls = map(x -> url_given_size(x, "square_medium"), current_page_illusts)
    return urls
end

function post_titles_in_page(current_page_illusts)
    titles = map(x -> post_title(current_page_illusts, x), 1:length(current_page_illusts))
    return titles
end

function page_urls_in_post(post_json, size="medium")
    """Get the number of pages and each of their urls in a multi-image post."""
    number_of_pages = post_json["page_count"]
    if number_of_pages > 1
        println("Page 1/$number_of_pages")
        list_of_pages = post_json["meta_pages"]
        page_urls = [url_given_size(list_of_pages[i], size)
                    for i in 1:number_of_pages]
        #for i in 1:number_of_pages
        #    append!(page_urls, url_given_size(list_of_pages[i], size))
        #end
    else
        page_urls = [url_given_size(post_json, size)]
    end

    return number_of_pages, page_urls
end

function change_url_to_full(url::String, png=false)
    """
    The 'large' resolution url isn't the largest. This uses changes the url to
    the highest resolution available
    """
    url = replace(url, r"_master\d+" => "")
    url = replace(url, r"c\/\d+x\d+_\d+_\w+\/img-master" => "img-original")

    # If it doesn't work, try changing to png
    if png
        url = replace(url, "jpg" => "png")
    end
    return url
end

function change_url_to_full(post_json, png=false)
    """
    The 'large' resolution url isn't the largest. This uses changes the url to
    the highest resolution available
    """
    if post_json
        url = url_given_size(post_json, "large")
    end
    url = replace(url, r"_master\d+" => "")
    url = replace(url, r"c\/\d+x\d+_\d+_\w+\/img-master" => "img-original")

    # If it doesn't work, try changing to png
    if png
        url = replace(url, "jpg" => "png")
    end
    return url
end

function process_user_url(url_or_id)
    if occursin("users", url_or_id)
        if occursin("\\", url_or_id)
            user_input = split(split_backslash_last(url_or_id), "\\")[end][2:end]
        else
            user_input = split_backslash_last(url_or_id)
        end
    else
        user_input = url_or_id
    end
    return user_input, "1"
end

function process_artwork_url(url_or_id)
    if occursin("artworks", url_or_id)
        user_input = split(split_backslash_last(url_or_id), "\\")[begin]
    elseif occursin("illust_id", url_or_id)
        user_input = match(r"&illust_id.*", url_or_id)#[0].split("=")[-1]
    else
        user_input = url_or_id
    end
    return user_input, "2"
end

end


text:
"""

  This function only sets the value in mixed layer without
  conserving quantity

"""
 function setMixedLayer!(;
    Ts   :: AbstractArray{Float64, 1},
    Ss   :: AbstractArray{Float64, 1},
    zs   :: AbstractArray{Float64, 1},
    Nz   :: Integer,
    T_ML :: Float64,
    S_ML :: Float64,
    h_ML :: Float64,
)
    FLDO = getFLDO(zs=zs, h_ML=h_ML, Nz=Nz)

    if FLDO > 1
        Ts[1:FLDO-1] .= T_ML
        Ss[1:FLDO-1] .= S_ML
    elseif FLDO == -1
        Ts[1:Nz] .= T_ML
        Ss[1:Nz] .= S_ML
    end

    return FLDO
end

function OC_setMixedLayer!(
    ocn  :: Ocean,
    i    :: Integer,
    j    :: Integer;
    T_ML :: Float64,
    S_ML :: Float64,
    h_ML :: Float64,
)

    ocn.h_ML[i, j] = h_ML
    ocn.T_ML[i, j] = T_ML
    ocn.S_ML[i, j] = S_ML
    ocn.FLDO[i, j] = setMixedLayer!(
        Ts   = ocn.cols.Ts[i, j],
        Ss   = ocn.cols.Ss[i, j],
        zs   = ocn.cols.zs[i, j],
        Nz   = ocn.Nz[i, j],
        T_ML = T_ML,
        S_ML = S_ML,
        h_ML = h_ML,
    )

end

function OC_setBuoyancy!(
    ocn  :: Ocean,
    i    :: Integer,
    j    :: Integer;
    bs   :: AbstractArray{Float64,1},
    b_ML :: Float64,
    h_ML :: Float64,
)

    ocn.bs[i, j, :] = bs
    OC_setMixedLayer!(ocn, i, j; b_ML=b_ML, h_ML=h_ML)

end

#=
function makeBlankOceanColumn(;zs::Array{Float64, 1})
    N     = length(zs) - 1
    bs    = zeros(Float64, N)
    K     = 0.0
    b_ML  = 0.0
    h_ML  = h_ML_min
    FLDO  = 1

    oc = OceanColumn(N=N, zs=zs, bs=bs, K=K, b_ML=b_ML, h_ML=h_ML, FLDO=FLDO)
    OC_updateFLDO!(oc)

    return oc
end

function makeSimpleOceanColumn(;
    zs      :: AbstractArray{Float64, 1},
    b_slope :: Float64 = 30.0 / 5000.0 * g * α,
    b_ML    :: Float64 = 1.0,
    h_ML    :: Float64 = h_ML_min,
    Δb      :: Float64 = 0.0,
    K       :: Float64 = 1e-5
)

    oc = makeBlankOceanColumn(zs=zs)

    bs = zeros(Float64, length(zs)-1)
    for i = 1:length(bs)
        z = (zs[i] + zs[i+1]) / 2.0
        if z > -h_ML
            bs[i] = b_ML
        else
            bs[i] = b_ML - Δb - b_slope * (-z - h_ML)
        end
    end

    OC_setBuoyancy!(oc, bs=bs, b_ML=b_ML, h_ML=h_ML)
    oc.K = K
    OC_updateFLDO!(oc)

    return oc

end
=#

max_stars_repo_name: meteorologytoday/IOM

max_issues_repo_name: meteorologytoday/IOM

max_forks_repo_name: meteorologytoday/IOM


text:
export stationary_initial_params_tilde, stationary_initial_params, dynamic_initial_params

"""
    stationary_initial_param_tilde(gas::Model{D, T}) where {D, T}

#TODO
"""
function stationary_initial_params_tilde(gas::Model{D, T}) where {D, T}
    biggest_lag = number_of_lags(gas)
    n_params = num_params(D)
    initial_params_tilde = Matrix{T}(undef, biggest_lag, n_params)
    for t in 1:biggest_lag, p in 1:n_params
        # ω/(1 - sum(B[p, p]))
        initial_params_tilde[t, p] = gas.ω[p]/(1 - sum(v[p, p] for (k, v) in gas.B))
    end
    return initial_params_tilde
end

"""
    stationary_initial_params(gas::Model{D, T}) where {D, T}

#TODO
"""
function stationary_initial_params(gas::Model{D, T}) where {D, T}
    biggest_lag = number_of_lags(gas)
    n_params = num_params(D)
    initial_params_tilde = Matrix{T}(undef, biggest_lag, n_params)
    initial_params = Matrix{T}(undef, biggest_lag, n_params)
    for t in 1:biggest_lag, p in 1:n_params
        # ω/(1 - sum(B[p, p]))
        initial_params_tilde[t, p] = gas.ω[p]/(1 - sum(v[p, p] for (k, v) in gas.B))
        unlink!(initial_params, D, initial_params_tilde, t)
    end
    return initial_params
end

"""
    dynamic_initial_params

    #TODO This is an heuristic! only seen in some phd thesis.
"""
function dynamic_initial_params(obs::Vector{T}, gas::Model{D, T}) where {D, T}
    # Take the biggest lag
    biggest_lag = ScoreDrivenModels.number_of_lags(gas)

    # Allocate memory
    initial_params = Matrix{T}(undef, biggest_lag, ScoreDrivenModels.num_params(D))
    obs_separated = Vector{Vector{T}}(undef, biggest_lag)

    # Loop to fit mle in every component of seasonality
    len = length(obs)
    for i in 1:biggest_lag
        # Indexes of a component of seasonality i.e. of january in a annual seasonality
        idx = collect(i:biggest_lag:len)
        # Fit MLE in the observations
        dist = fit_mle(D, obs[idx])
        # Adequate to the ScoreDrivenModels standard
        # In Distributions Normal is \mu and \sigma
        # In ScoreDrivenModels Normal is \mu and \sigma^2
        initial_params[i, :] = permutedims([ScoreDrivenModels.params_sdm(dist)...])
    end

    return initial_params
end


text:
@testset "RandomObs (and overlaps with BalancedObs)" begin
    for TIter in (RandomObs, BalancedObs)
        @test TIter <: LearnBase.DataIterator
        @test TIter <: LearnBase.ObsIterator
        @test TIter <: LearnBase.AbstractDataIterator
        @test TIter <: LearnBase.AbstractObsIterator
    end
    @test_reference "references/RandomObs1.txt" RandomObs(X)
    @test_reference "references/RandomObs2.txt" RandomObs(X, 10)
    @test_reference "references/BalancedObs1.txt" BalancedObs(X)
    @test_reference "references/BalancedObs2.txt" BalancedObs(X, 10)
    if  Int == Int64
        @test_reference "references/RandomObs3.txt" [RandomObs(X)]
        @test_reference "references/BalancedObs3.txt" [BalancedObs(X)]
    end

    for TIter in (RandomObs, BalancedObs)
        @testset "constructor for $TIter" begin
            for var in (vars..., tuples..., Xs, ys)
                A = @inferred TIter(var)
                @test_throws MethodError collect(A)
                @test typeof(A) <: TIter
                @test @inferred(nobs(A)) == 150
                @test_throws MethodError length(A)
                @test typeof(@inferred(first(A))) == typeof(datasubset(var, 1))
                A = @inferred TIter(var, 100)
                @test typeof(A) <: TIter
                @test @inferred(nobs(A)) == 150
                @test @inferred(length(A)) == 100
                @test typeof(@inferred(first(A))) == typeof(datasubset(var, 1))
                @test eltype(@inferred(collect(A))) == typeof(datasubset(var, 1))
                B = TIter(var, count = 100)
                @test typeof(B) == typeof(A)
                @test A.count == B.count
            end
            for var in (X,XX,XXX)
                A = @inferred TIter(var, 100)
                @test size(getobs(first(A))) == size(getobs(datasubset(var, 1)))
            end
        end
    end

    for TIter in (RandomObs, BalancedObs)
        @testset "DataSubset: $TIter" begin
            for raw in (vars..., tuples..., Xs, ys)
                var = DataSubset(raw, 1:90)
                A = @inferred TIter(var)
                @test_throws MethodError collect(A)
                @test typeof(A) <: TIter
                @test @inferred(nobs(A)) == 90
                @test_throws MethodError length(A)
                @test_throws MethodError getobs(A)
                @test_throws MethodError getobs(A, 1)
                @test typeof(@inferred(first(A))) == typeof(datasubset(var, 1))
                A = @inferred TIter(var, 100)
                @test typeof(A) <: TIter
                @test @inferred(nobs(A)) == 90
                @test @inferred(length(A)) == 100
                @test_throws MethodError getobs(A, 1)
                @test typeof(@inferred(first(A))) == typeof(datasubset(var, 1))
                @test eltype(@inferred(collect(A))) == typeof(datasubset(var, 1))
                @test typeof(@inferred(getobs(A))) == typeof(getobs.(collect(A)))
            end
        end
    end

    for TIter in (RandomObs, BalancedObs)
        @testset "various obsdim: $TIter" begin
            A = @inferred TIter(X', ObsDim.First())
            B = TIter(X', obsdim = 1)
            @test A.data == B.data
            @test A.obsdim == B.obsdim
            @test A.count == B.count
            A = @inferred TIter(X', 10, ObsDim.First())
            B = TIter(X', 10, obsdim = 1)
            C = TIter(X', count = 10, obsdim = 1)
            @test A.data == B.data == C.data
            @test A.obsdim == B.obsdim == C.obsdim
            @test A.count == B.count == C.count
            @test typeof(@inferred(first(A))) == typeof(datasubset(X', 1, ObsDim.First()))
            @test eltype(@inferred(collect(A))) == typeof(datasubset(X', 1, ObsDim.First()))
            @test size(getobs(first(A))) == size(getobs(datasubset(X', 1, ObsDim.First())))
            A = @inferred TIter((X',X'), 10, ObsDim.First())
            @test typeof(@inferred(first(A))) == typeof(datasubset((X',X'), 1, ObsDim.First()))
            @test eltype(@inferred(collect(A))) == typeof(datasubset((X',X'), 1, ObsDim.First()))
            A = @inferred TIter((X',X), 10, (ObsDim.First(),ObsDim.Last()))
            @test typeof(@inferred(first(A))) == typeof(datasubset((X',X), 1, (ObsDim.First(),ObsDim.Last())))
            @test eltype(@inferred(collect(A))) == typeof(datasubset((X',X), 1, (ObsDim.First(),ObsDim.Last())))
            A = @inferred TIter(datasubset(X', ObsDim.First()))
        end
    end

    for TIter in (RandomObs, BalancedObs)
        @testset "custom types: TIter" begin
            @test_throws MethodError TIter(EmptyType())
            @test_throws MethodError TIter(EmptyType(), 10)
            @test_throws MethodError TIter(EmptyType(), obsdim=1)
            @test_throws MethodError TIter(EmptyType(), 10, obsdim=1)
            A = @inferred TIter(CustomType())
            @test_throws MethodError collect(A)
            @test typeof(A) <: TIter
            @test @inferred(nobs(A)) == 100
            @test_throws MethodError length(A)
            @test typeof(@inferred(first(A))) == typeof(datasubset(CustomType(), 1))
            A = @inferred TIter(CustomType(), 500)
            @test typeof(A) <: TIter
            @test @inferred(nobs(A)) == 100
            @test @inferred(length(A)) == 500
            @test typeof(@inferred(first(A))) == typeof(datasubset(CustomType(), 1))
            @test eltype(@inferred(collect(A))) == typeof(datasubset(CustomType(), 1))
            @test all(getobs.(collect(A)) .<= 100)
        end
    end
end

@testset "BalancedObs" begin
    @testset "check balancing" begin
        A = @inferred BalancedObs((1:10,[:a,:b,:a,:a,:a,:a,:a,:a,:b,:a]), 100)
        @test A.labelmap == Dict(:a => [1,3,4,5,6,7,8,10], :b => [2,9])
        val = @inferred(map(getobs, A))
        freq = labelfreq(getindex.(val, 2))
        @test freq[:a] > 35
        @test freq[:b] > 35

        A = @inferred BalancedObs(i->i>1, 1:10, 100)
        @test A.labelmap == Dict(false => [1], true => [2,3,4,5,6,7,8,9,10])
        val = @inferred(map(getobs, A))
        @test sum(val .> 1) > 35
    end
end

@testset "RandomBatches" begin
    @test RandomBatches <: LearnBase.DataIterator
    @test RandomBatches <: LearnBase.BatchIterator
    @test RandomBatches <: LearnBase.AbstractDataIterator
    @test RandomBatches <: LearnBase.AbstractBatchIterator
    @test_reference "references/RandomBatches1.txt" RandomBatches(X)
    @test_reference "references/RandomBatches2.txt" RandomBatches(X,10,10)
    if Int == Int64
        @test_reference "references/RandomBatches3.txt" [RandomBatches(X)]
    end

    @testset "constructor" begin
        A = @inferred RandomBatches(rand(2,5), 10)
        @test typeof(A) <: RandomBatches
        @test @inferred(nobs(A)) == 5
        @test @inferred(batchsize(A)) == 10
        @test typeof(@inferred(first(A))) == typeof(datasubset(rand(2,5), rand(1:5, 10)))

        for var in (vars..., tuples..., Xs, ys)
            A = @inferred RandomBatches(var)
            @test_throws MethodError collect(A)
            @test typeof(A) <: RandomBatches
            @test @inferred(nobs(A)) == 150
            @test @inferred(batchsize(A)) == 30
            @test_throws MethodError length(A)
            @test typeof(@inferred(first(A))) == typeof(datasubset(var, collect(1:30)))
            A = @inferred RandomBatches(var, 50)
            @test_throws MethodError collect(A)
            @test typeof(A) <: RandomBatches
            @test @inferred(nobs(A)) == 150
            @test @inferred(batchsize(A)) == 50
            @test_throws MethodError length(A)
            @test typeof(@inferred(first(A))) == typeof(datasubset(var, collect(1:50)))
            B = RandomBatches(var, size = 50)
            @test typeof(A) == typeof(B)
            @test A.data == B.data
            @test A.count == B.count
            @test A.size == B.size
            @test A.obsdim == B.obsdim
            A = @inferred RandomBatches(var, 40, 100)
            @test typeof(A) <: RandomBatches
            @test @inferred(nobs(A)) == 150
            @test @inferred(batchsize(A)) == 40
            @test @inferred(length(A)) == 100
            @test typeof(@inferred(first(A))) == typeof(datasubset(var, collect(1:40)))
            @test eltype(@inferred(collect(A))) == typeof(datasubset(var, collect(1:40)))
            B = RandomBatches(var, size = 40, count = 100)
            @test typeof(A) == typeof(B)
            @test A.data == B.data
            @test A.count == B.count
            @test A.size == B.size
            @test A.obsdim == B.obsdim
        end
        for var in (X,XX,XXX)
            A = @inferred RandomBatches(var)
            @test size(getobs(first(A))) == size(getobs(datasubset(var, collect(1:30))))
            A = @inferred RandomBatches(var, 50)
            @test size(getobs(first(A))) == size(getobs(datasubset(var, collect(1:50))))
            A = @inferred RandomBatches(var, 40, 100)
            @test size(getobs(first(A))) == size(getobs(datasubset(var, collect(1:40))))
        end
    end

    @testset "DataSubset" begin
        for raw in (vars..., tuples..., Xs, ys)
            var = DataSubset(raw)
            A = @inferred RandomBatches(var)
            @test_throws MethodError collect(A)
            @test typeof(A) <: RandomBatches
            @test @inferred(nobs(A)) == 150
            @test @inferred(batchsize(A)) == 30
            @test_throws MethodError length(A)
            @test typeof(@inferred(first(A))) == typeof(datasubset(var, collect(1:30)))
            A = @inferred RandomBatches(var, 50)
            @test_throws MethodError collect(A)
            @test typeof(A) <: RandomBatches
            @test @inferred(nobs(A)) == 150
            @test @inferred(batchsize(A)) == 50
            @test_throws MethodError length(A)
            @test typeof(@inferred(first(A))) == typeof(datasubset(var, collect(1:50)))
            A = @inferred RandomBatches(var, 40, 100)
            @test typeof(A) <: RandomBatches
            @test @inferred(nobs(A)) == 150
            @test @inferred(batchsize(A)) == 40
            @test @inferred(length(A)) == 100
            @test typeof(@inferred(first(A))) == typeof(datasubset(var, collect(1:40)))
            @test eltype(@inferred(collect(A))) == typeof(datasubset(var, collect(1:40)))
        end
    end

    @testset "various obsdim" begin
        A = @inferred RandomBatches(X', 10, ObsDim.First())
        B = RandomBatches(X', 10, obsdim = 1)
        @test A.data == B.data
        @test A.obsdim == B.obsdim
        @test A.count == B.count
        @test A.size == B.size == batchsize(A)
        A = @inferred RandomBatches(X', 10, 100, ObsDim.First())
        B = RandomBatches(X', 10, 100, obsdim = 1)
        C = RandomBatches(X', size = 10, count = 100, obsdim = 1)
        @test A.data == B.data == C.data
        @test A.obsdim == B.obsdim == C.obsdim
        @test A.count == B.count == C.count
        @test A.size == B.size == C.size == batchsize(A)
        @test typeof(@inferred(first(A))) == typeof(datasubset(X', collect(1:10), ObsDim.First()))
        @test eltype(@inferred(collect(A))) == typeof(datasubset(X', collect(1:10), ObsDim.First()))
        @test size(getobs(first(A))) == size(getobs(datasubset(X', collect(1:10), ObsDim.First())))
        A = @inferred RandomBatches(DataSubset(X', ObsDim.First()), 10, 100)
        @test size(getobs(first(A))) == size(getobs(datasubset(X', collect(1:10), ObsDim.First())))
    end

    @testset "custom types" begin
        @test_throws MethodError RandomBatches(EmptyType())
        @test_throws MethodError RandomBatches(EmptyType(), 10)
        @test_throws MethodError RandomBatches(EmptyType(), 10, 100)
        @test_throws MethodError RandomBatches(EmptyType(), obsdim=1)
        @test_throws MethodError RandomBatches(EmptyType(), 10, obsdim=1)
        A = @inferred RandomBatches(CustomType(), 20)
        @test_throws MethodError collect(A)
        @test typeof(A) <: RandomBatches
        @test @inferred(nobs(A)) == 100
        @test @inferred(batchsize(A)) == 20
        @test_throws MethodError length(A)
        @test typeof(@inferred(first(A))) == typeof(datasubset(CustomType(), collect(1:20)))
        A = @inferred RandomBatches(CustomType(), 20, 500)
        @test typeof(A) <: RandomBatches
        @test @inferred(nobs(A)) == 100
        @test @inferred(length(A)) == 500
        @test @inferred(batchsize(A)) == 20
        @test typeof(@inferred(first(A))) == typeof(datasubset(CustomType(), collect(1:20)))
        @test eltype(@inferred(collect(A))) == typeof(datasubset(CustomType(), collect(1:20)))
    end
end

@testset "BufferGetObs" begin
    @testset "ObsView" begin
        A = BufferGetObs(ObsView(X))
        @test_reference "references/BufferGetObs1.txt" A
        if Int == Int64
            @test_reference "references/BufferGetObs2.txt" [A]
        end

        @test size(A.buffer) == (4,)
        @test typeof(A.buffer) <: Array{Float64,1}

        for var in (X,Y,XX,XXX,(X,Y),(XX,X),(XXX,XX,X))
            A = @inferred BufferGetObs(ObsView(var))
            @test typeof(A.buffer) == typeof(getobs(var,1))
            for (i,a) in enumerate(A)
                @test a === A.buffer
                @test a == getobs(var,i)
            end
        end

        # preallocating buffer
        for var in (X,Y,XX,XXX)
            buffer = similar(getobs(var,1))
            A = @inferred BufferGetObs(ObsView(var), buffer)
            for (i,a) in enumerate(A)
                @test a === A.buffer
                @test a === buffer
                @test a == getobs(var,i)
            end
        end
        # preallocating buffer for tuple
        for var in ((X,X),(XX,X),(X,Y))
            buffer = (similar(getobs(var[1],1)),similar(getobs(var[2],1)))
            A = @inferred BufferGetObs(ObsView(var), buffer)
            for (i,a) in enumerate(A)
                @test a === A.buffer
                @test a === buffer
                @test a == getobs(var,i)
            end
        end

        for var in (vars..., tuples..., Xs, ys, CustomType())
            A = @inferred BufferGetObs(ObsView(var))
            @test typeof(A.buffer) == typeof(getobs(var,1))
            for (i,a) in enumerate(A)
                @test a == getobs(var,i)
            end
        end
    end

    @testset "BatchView" begin
        A = BufferGetObs(BatchView(X, 15))
        @test size(A.buffer) == (4,15)
        @test typeof(A.buffer) <: Array{Float64,2}

        for var in (X,Y,XX,XXX,(X,Y),(XX,X),(XXX,XX,X))
            A = @inferred BufferGetObs(BatchView(var,15))
            @test typeof(A.buffer) == typeof(getobs(var,collect(1:15)))
            for (i,a) in enumerate(A)
                @test a === A.buffer
                @test a == BatchView(var,15)[i]
            end
        end
        # preallocating buffer
        for var in (X,Y,XX,XXX)
            buffer = similar(getobs(var,1:15))
            A = @inferred BufferGetObs(BatchView(var,15), buffer)
            for (i,a) in enumerate(A)
                @test a === A.buffer
                @test a === buffer
                @test a == BatchView(var,15)[i]
            end
        end
        # preallocating buffer for tuple
        for var in ((X,X),(XX,X),(X,Y))
            buffer = (similar(getobs(var[1],1:15)),similar(getobs(var[2],1:15)))
            A = @inferred BufferGetObs(BatchView(var, 15), buffer)
            for (i,a) in enumerate(A)
                @test a === A.buffer
                @test a === buffer
                @test a == BatchView(var,15)[i]
            end
        end

        for var in (vars..., tuples..., Xs, ys, CustomType())
            A = @inferred BufferGetObs(BatchView(var,10))
            @test typeof(A.buffer) == typeof(getobs(var,1:10))
            for (i,a) in enumerate(A)
                @test a == getobs(BatchView(var,10)[i])
            end
        end
    end

    @testset "RandomObs" begin
        A = BufferGetObs(RandomObs(X, 10))
        @test size(A.buffer) == (4,)
        @test typeof(A.buffer) <: Array{Float64,1}
        @test length(A) == 10
        @test Base.IteratorSize(A) == Base.HasLength()
        A = BufferGetObs(RandomObs(X))
        @test size(A.buffer) == (4,)
        @test typeof(A.buffer) <: Array{Float64,1}
        @test_throws MethodError length(A)
        @test Base.IteratorSize(A) == Base.IsInfinite()

        for var in (X,Y,XX,XXX,(X,Y),(XX,X),(XXX,XX,X))
            A = @inferred BufferGetObs(RandomObs(var,10))
            @test typeof(A.buffer) == typeof(getobs(var,1))
            for (i,a) in enumerate(A)
                @test a === A.buffer
                @test i <= 10
            end
        end
        for var in (vars..., tuples..., Xs, ys, CustomType())
            A = @inferred BufferGetObs(RandomObs(var, 10))
            @test typeof(A.buffer) == typeof(getobs(var,1))
            for (i,a) in enumerate(A)
                @test i <= 10
            end
        end
    end

    @testset "RandomBatches" begin
        A = BufferGetObs(RandomBatches(X, 15, 10))
        @test size(A.buffer) == (4,15)
        @test typeof(A.buffer) <: Array{Float64,2}
        @test length(A) == 10
        @test Base.IteratorSize(A) == Base.HasLength()
        A = BufferGetObs(RandomBatches(X, 15))
        @test size(A.buffer) == (4,15)
        @test typeof(A.buffer) <: Array{Float64,2}
        @test_throws MethodError length(A)
        @test Base.IteratorSize(A) == Base.IsInfinite()

        for var in (X,Y,XX,XXX,(X,Y),(XX,X),(XXX,XX,X))
            A = @inferred BufferGetObs(RandomBatches(var,15,10))
            @test typeof(A.buffer) == typeof(getobs(var,collect(1:15)))
            for (i,a) in enumerate(A)
                @test a === A.buffer
                @test i <= 10
            end
        end
        for var in (vars..., tuples..., Xs, ys, CustomType())
            A = @inferred BufferGetObs(RandomBatches(var,10,10))
            @test typeof(A.buffer) == typeof(getobs(var,1:10))
            for (i,a) in enumerate(A)
                @test i <= 10
            end
        end
    end
end

@testset "eachobs" begin
    A = @inferred eachobs(X)
    @inferred first(A)
    @test typeof(A) <: BufferGetObs
    @test typeof(A.iter) <: ObsView
    @test A.iter.obsdim == ObsDim.Last()
    @test A.iter.data == X
    A = @inferred eachobs(X, ObsDim.First())
    @inferred first(A)
    @test typeof(A) <: BufferGetObs
    @test typeof(A.iter) <: ObsView
    @test A.iter.obsdim == ObsDim.First()
    @test A.iter.data == X
    A = eachobs(X, obsdim = 1)
    @inferred first(A)
    @test typeof(A) <: BufferGetObs
    @test typeof(A.iter) <: ObsView
    @test A.iter.obsdim == ObsDim.First()
    @test A.iter.data == X
end

@testset "eachbatch" begin
    A = @inferred eachbatch(X)
    @inferred first(A)
    @test typeof(A) <: BufferGetObs
    @test typeof(A.iter) <: BatchView
    @test A.iter.obsdim == ObsDim.Last()
    @test A.iter.data == X
    @test size(A.buffer) == (4,30)
    @test length(A) == 5
    A = @inferred eachbatch(X', ObsDim.First())
    @inferred first(A)
    @test typeof(A) <: BufferGetObs
    @test typeof(A.iter) <: BatchView
    @test A.iter.obsdim == ObsDim.First()
    @test A.iter.data == X'
    @test size(A.buffer) == (30,4)
    @test length(A) == 5
    A = eachbatch(X', obsdim=1)
    @inferred first(A)
    @test typeof(A) <: BufferGetObs
    @test typeof(A.iter) <: BatchView
    @test A.iter.obsdim == ObsDim.First()
    @test A.iter.data == X'
    @test size(A.buffer) == (30,4)
    @test length(A) == 5
    A = @inferred eachbatch(X, 10)
    @test_throws ArgumentError eachbatch(X, size = 10, maxsize = 11)
    @test length(A) == length(eachbatch(X, maxsize = 11))
    @inferred first(A)
    @test typeof(A) <: BufferGetObs
    @test typeof(A.iter) <: BatchView
    @test A.iter.obsdim == ObsDim.Last()
    @test A.iter.data == X
    @test size(A.buffer) == (4,10)
    @test length(A) == 15
    A = @inferred eachbatch(X', 10, ObsDim.First())
    @inferred first(A)
    @test typeof(A) <: BufferGetObs
    @test typeof(A.iter) <: BatchView
    @test A.iter.obsdim == ObsDim.First()
    @test A.iter.data == X'
    @test size(A.buffer) == (10,4)
    @test length(A) == 15
    A = eachbatch(X', count = 15, obsdim = 1)
    @inferred first(A)
    @test typeof(A) <: BufferGetObs
    @test typeof(A.iter) <: BatchView
    @test A.iter.obsdim == ObsDim.First()
    @test A.iter.data == X'
    @test size(A.buffer) == (10,4)
    @test length(A) == 15
end


text:
#   Unit examples of the FourierAnalysis Package for julia language
#   v 0.0.1 - last update 24th of September 2019
#
#   MIT License
#   Copyright (c) 2019, Marco Congedo, CNRS, Grenobe, France:
#   https://sites.google.com/site/marcocongedo/home

# ? CONTENTS :
#   This example show how to compute the analytic signal via
#   Hilbert transform and how to extract inforation from them.

#   ~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~  #
#                                                                             #
#   ~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~~¤~  #

using FourierAnalysis, FFTW, LinearAlgebra, Statistics, Plots

# add module for reading the two EEG text files to be used ater
push!(LOAD_PATH, @__DIR__)
using IOtxt

# Get EEG file names with complete path (they have extension .txt)
S=getFilesInDir(@__DIR__; ext=(".txt",))

# read the two EEG data files and put them in a Matrix object
X1=readEEG(S[1])
X2=readEEG(S[2])

######################################################################
# Test Hilbert Transform Analytic Signal

using FourierAnalysis, FFTW, LinearAlgebra, Statistics, Plots, DSP
t=128; lab=["x", "real(y)", "imag(y)"]

# Analytic signal of one vector
x=sinusoidal(10, 2, 128, t, π/2; DC=10) # cosine
y=analyticsignal(x)
# note that real(y) is x without the DC level, i.e., x=real(y)+DC
plot([x, real(y), imag(y)]; labels=lab)

# make a check
s=sinusoidal(10, 2, 128, t, 0) # sine
norm(s-imag(y)) # should be zero

# Analytic Signal by DSP.jl
y2=hilbert(x)
norm(s-imag(y2)) # should be zero
# DSP.jl does not remove the DC level
# thus x=real(y2) in this case
plot([x, real(y2), imag(y2)]; labels=lab)

# Analytic signal of multiple vectors
x=hcat(x, sinusoidal(10, 3, 128, t, π/2; DC=10))
y=analyticsignal(x)

# welch-like analytic signal of one vector
# (note edge effects)
x=sinusoidal(10, 2, 128, t*4, π/2; DC=0)
y=analyticsignal(x, t)
plot([x, real(y), imag(y)]; labels=lab)

# Welch-like analytic signal of multiple vectors
x=hcat(x, sinusoidal(10, 3, 128, t*4, π/2; DC=0))
y=analyticsignal(x, t)

t=256
lim=512
Y=analyticsignal(X1[1:lim, :], lim)
plot([X1[1:lim, 1], real(Y[:, 1]), imag(Y[:, 1])])

# compute Welch-like Analytis Signal altogether for X1 and X2
𝒀 = analyticsignal([X1[1:lim, :], X2[1:lim, :]], lim)

## test amplitude, phase and polar functions

x=sinusoidal(10, 2, 128, t*4, 0).*sinusoidal(10, 1, 128, t*4, 0)

# amplitude and phase of a vector using analytic signal standard method
y=analyticsignal(x)
a=amplitude(y)
ϕ=phase(y, func=x->(x+π)/2π*50)
plot(x; labels="signal")
plot!(a; labels="amplitude")
plot!(ϕ; labels="phase")

# see what happen if `x` contains energy in frequencies below sr/wl Hz
# (see documentation of `analyticSignal` function)
y=analyticsignal(x, 64)
a=amplitude(y)
ϕ=phase(y, func=x->(x+π)/2π*50)
plot(x; labels="signal")
plot!(a; labels="amplitude")
plot!(ϕ; labels="phase")

# get Amplitude from analytic Signal of a data matrix holding multiple time-series
A=amplitude(Y)
plot(A[:, 1:2])

# get Phase from analytic Signal
𝛷=phase(Y)
plot(𝛷[:, 1:1])

# get Phase from analytic Signal and transform it in [-1, 1]
𝛷=phase(Y, func=x->(x+π)/2π)
plot(𝛷[:, 1:1])

# get Phase from analytic Signal and transform it to the sine of the phase
𝛷=phase(Y, func=x->sin(x))
plot(𝛷[:, 1:1])

# get Amplitude and phase from analytic Signal
A, 𝛷=polar(Y)


text:
export cg1

"""
`cg1(j1, m1, j2, m2, j3, m3, T=Float64)` : A reference implementation of
Clebsch-Gordon coefficients based on

https://hal.inria.fr/hal-01851097/document
Equation (4-6)

This heavily uses BigInt and BigFloat and should therefore not be employed
for performance critical tasks.
"""
function cg1(j1, m1, j2, m2, j3, m3, T=Float64)
   if (m3 != m1 + m2) || !(abs(j1-j2) <= j3 <= j1 + j2)
      return zero(T)
   end

   N = (2*j3+1) *
       factorial(big(j1+m1)) * factorial(big(j1-m1)) *
       factorial(big(j2+m2)) * factorial(big(j2-m2)) *
       factorial(big(j3+m3)) * factorial(big(j3-m3)) /
       factorial(big( j1+j2-j3)) /
       factorial(big( j1-j2+j3)) /
       factorial(big(-j1+j2+j3)) /
       factorial(big(j1+j2+j3+1))

   G = big(0)
   # 0 ≦ k ≦ j1+j2-j3
   # 0 ≤ j1-m1-k ≤ j1-j2+j3   <=>   j2-j3-m1 ≤ k ≤ j1-m1
   # 0 ≤ j2+m2-k ≤ -j1+j2+j3  <=>   j1-j3+m2 ≤ k ≤ j2+m2
   lb = (0, j2-j3-m1, j1-j3+m2)
   ub = (j1+j2-j3, j1-m1, j2+m2)
   for k in maximum(lb):minimum(ub)
      bk = big(k)
      G += (-1)^k *
           binomial(big( j1+j2-j3), big(k)) *
           binomial(big( j1-j2+j3), big(j1-m1-k)) *
           binomial(big(-j1+j2+j3), big(j2+m2-k))
   end

   return T(sqrt(N) * G)
end


text:
# This file was generated by the Julia Swagger Code Generator
# Do not modify this file directly. Modify the swagger specification instead.

mutable struct TagDetails <: SwaggerModel
    id::Any # spec type: Union{ Nothing, String } # spec name: id
    tagName::Any # spec type: Union{ Nothing, String } # spec name: tagName
    count::Any # spec type: Union{ Nothing, TagCount } # spec name: count
    values::Any # spec type: Union{ Nothing, Vector{TagValue} } # spec name: values

    function TagDetails(;id=nothing, tagName=nothing, count=nothing, values=nothing)
        o = new()
        validate_property(TagDetails, Symbol("id"), id)
        setfield!(o, Symbol("id"), id)
        validate_property(TagDetails, Symbol("tagName"), tagName)
        setfield!(o, Symbol("tagName"), tagName)
        validate_property(TagDetails, Symbol("count"), count)
        setfield!(o, Symbol("count"), count)
        validate_property(TagDetails, Symbol("values"), values)
        setfield!(o, Symbol("values"), values)
        o
    end
end # type TagDetails

const _property_map_TagDetails = Dict{Symbol,Symbol}(Symbol("id")=>Symbol("id"), Symbol("tagName")=>Symbol("tagName"), Symbol("count")=>Symbol("count"), Symbol("values")=>Symbol("values"))
const _property_types_TagDetails = Dict{Symbol,String}(Symbol("id")=>"String", Symbol("tagName")=>"String", Symbol("count")=>"TagCount", Symbol("values")=>"Vector{TagValue}")
Base.propertynames(::Type{ TagDetails }) = collect(keys(_property_map_TagDetails))
Swagger.property_type(::Type{ TagDetails }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_TagDetails[name]))}
Swagger.field_name(::Type{ TagDetails }, property_name::Symbol) =  _property_map_TagDetails[property_name]

function check_required(o::TagDetails)
    true
end

function validate_property(::Type{ TagDetails }, name::Symbol, val)
end