text:
ENV["COLUMNS"]=72
import Base.GC.gc
using DSP.Periodograms: stft, nextfastfft
using DSP.Windows: hanning
using HDF5
using Base.Iterators: flatten, cycle, take, partition
using Statistics: mean, std
using Random: shuffle!
using Knet
using Knet: Data, KnetArray
using LinearAlgebra: diag

function epsnormalise(x::AbstractArray{T}; dims=1, ϵ=eps(T)) where T
    μ = mean(x, dims = dims)
    σ = std(x, dims = dims, mean = μ, corrected=false)
    σ = map((x) -> abs(x) < ϵ ? ϵ : x, σ)
    return (x .- μ) ./ σ
end

function make_batch(data, labels, idxs, minibatch_size, cat_mapping; divider = 1, n = 512, noverlap = div(n, 2))
    nfft = nextfastfft(n)
    sz = div(div(size(data, 1), noverlap), divider) - 1
    X_batch = Array{Float32}(undef, div(nfft, 2) + 1, sz, 1, length(idxs))
    Y_batch = Array{UInt8}(undef, length(idxs))
    win = hanning(nfft)
    data_ind = 1:div(size(data, 1), divider)
    for i in 1:length(idxs)
        st = abs.(stft(Float32.(data[data_ind,idxs[i]][:]./typemax(Int16)), n, noverlap, window=win, nfft=nfft))
        X_batch[:, :, 1, i] = epsnormalise(20*log10.(max.(st, eps(Float32))))
        Y_batch[i] = cat_mapping[labels[idxs[i]][1]]
    end
    return minibatch(convert(KnetArray, X_batch), Y_batch, minibatch_size)
end

function confusion(model, data, no_categories, no_labels=no_categories)
    conf = zeros(no_categories, no_labels)
    sz = 0
    for mbatch in data
        pred = convert(Array, model(mbatch[1]))
        for ind = 1:length(mbatch[2])
            conf[findmax(pred[:, ind])[2], mbatch[2][ind]] += 1
        end
    end
    return conf
end

function squish_confusion(conf_ex, cat_mapping)
    vals = sort(collect(values(cat_mapping)))
    no_categories = length(unique(vals))
    conf = zeros(no_categories, no_categories)
    for i = 1:length(vals)
        conf[:, vals[i]] += conf_ex[:, i]
    end
    conf
end

function get_idxs(size, batch_size, cat_mapping, start_idx=1)
    idxs = Vector{Vector{Int}}()
    no_categories = length(unique(values(cat_mapping)))
    no_labels = length(cat_mapping)
    cat_lengths = zeros(no_categories)
    for val in values(cat_mapping)
        cat_lengths[val] += 1
    end
    cat_idxs = zeros(no_categories)

    batch_i = 1
    lab_i = start_idx
    push!(idxs, Vector{UInt}())
    while batch_i <= size/batch_size
        for i = 1:length(cat_lengths)
            push!(idxs[batch_i], lab_i + sum(cat_lengths[1:(i-1)]) + cat_idxs[i])
            cat_idxs[i] = mod(cat_idxs[i] + 1, cat_lengths[i])
            if length(idxs[batch_i]) == batch_size
                batch_i += 1
                if batch_i > size/batch_size
                    break
                else
                    push!(idxs, Vector{UInt}())
                end
            end
        end
        lab_i += no_labels
    end
    idxs
end

cat_mapping =
Dict{Int8, UInt8}(0  => 1,  # WAV
                  11 => 2,  # MP3 320
                  12 => 3,  # MP3 192
                  13 => 4,  # MP3 128
                  21 => 5,  # AAC 320
                  22 => 6,  # AAC 192
                  23 => 7,  # AAC 128
                  31 => 8,  # OGG 320
                  32 => 9,  # OGG 192
                  33 => 10, # OGG 128
                  41 => 11, # WMA 320
                  42 => 12, # WMA 192
                  43 => 13, # WMA 128
                  51 => 14, # AC3 320
                  52 => 15, # AC3 192
                  53 => 16, # AC3 128
                 )

minibatch_size = 2
train_batch_size = 16
test_batch_size = 512
train_size = 5120
test_size = 2560

all_mapping =
Dict{Int8, UInt8}([(key, i) for (i, key) in enumerate(sort(collect(keys(cat_mapping))))])

no_categories = length(unique(values(cat_mapping)))
no_labels = length(cat_mapping)

train_batch_idxs = get_idxs(train_size, train_batch_size, cat_mapping)
test_batch_idxs = get_idxs(test_size, test_batch_size, cat_mapping, ceil(Int, train_batch_idxs[end][end]/no_labels)*no_labels + 1)

struct Conv; w; b; f; p; end
(c::Conv)(x) = c.f.(pool(conv4(c.w, dropout(x,c.p)) .+ c.b))
Conv(w1::Int,w2::Int,cx::Int,cy::Int,f=relu;pdrop=0) = Conv(param(w1,w2,cx,cy), param0(1,1,cy,1), f, pdrop)

struct Dense; w; b; f; p; end
(d::Dense)(x) = d.f.(d.w * mat(dropout(x,d.p)) .+ d.b) # mat reshapes 4-D tensor to 2-D matrix so we can use matmul
Dense(i::Int,o::Int,f=relu;pdrop=0) = Dense(param(o,i), param0(o), f, pdrop)

struct Chain; layers; Chain(layers...) = new(layers); end
(c::Chain)(x) = (for l in c.layers; x = l(x); end; x)
(c::Chain)(x,y) = nll(c(x),y)
(c::Chain)(d::Data) = mean(c(x,y) for (x,y) in d)

function train!(epochs, load=true, filename="model.jld2", train=true, test=true)
    gc()
    Knet.gc()

    model, test_confs = if load
        lf = Knet.load(filename)
        lf["model"], lf["test_confs"]
    else
        Chain(Conv(3, 3, 1, 16),
              Conv(3, 3, 16, 16),
              Conv(3, 3, 16, 16),
              Dense(102240, 300, pdrop=0.3),
              Dense(300, 300, pdrop=0.2),
              Dense(300, no_categories, identity, pdrop=0.2)),
        Array{Array{Float32, 2}, 1}()
    end
    best_acc = length(test_confs) > 0 ? sum(diag(test_confs[end]))/test_size : 0

    h5open("./dataset.h5") do file
        data = file["data"]
        labels = file["labels"]

        for epoch in 1:epochs
            @info "Epoch: $epoch"
            shuffle!(train_batch_idxs)
            if train
                for idx in progress(train_batch_idxs)
                    let dtrn = make_batch(data, labels, idx, minibatch_size, cat_mapping)
                        for a in adam(model, dtrn); end
                    end
                end
            end

            if test
                test_conf = zeros(no_categories, no_labels)
                for idx in progress(test_batch_idxs)
                    let dtst = make_batch(data, labels, idx, minibatch_size, all_mapping)
                        test_conf .+= confusion(model, dtst, no_categories, no_labels)
                    end
                end
                push!(test_confs, test_conf)

                acc = sum(diag(squish_confusion(test_conf, cat_mapping)))/test_size

                println("Confusion Matrix:")
                display(test_conf)
                println("Accuracy: $acc")
            else
                acc = best_acc + 1
            end
            if train
                if acc > best_acc
                    best_acc = acc
                    Knet.save(filename, "model", model, "test_confs", test_confs)
                end
            end
        end
    end

    return model, test_confs
end

max_stars_repo_name: Ashymad/IOLA

max_issues_repo_name: Ashymad/IOLA

max_forks_repo_name: Ashymad/IOLA


text:
function EuropeanCallTrinomial(K, T, S, sig, r, div, N, dx)

  dt = T/N
  nu = r - div - 0.5 * sig^2
  edx = exp(dx)
  pu = 0.5 * ((sig^2 * dt + nu^2 * dt^2) / dx^2 + nu * dt / dx)
  pm = 1 - (sig^2 * dt + nu^2 * dt^2) / dx^2
  pd = 0.5 * ((sig^2 * dt + nu^2 * dt^2) / dx^2 - nu * dt / dx)
  disc = exp(-r * dt)

  ## Initialise asset prices at maturity

  St = zeros(Float64, 2N+1)

  ## creating a shifted ranged for addressing the Array in the usual i,j notation
  J = sub(1:(2N+1), (1+(N+1):((2N+1)+(N+1))))

  ## Asset price at node N,-N
  St[J[-N]] = S * exp(-N * dx)

  ## Asset prices from bottom -N, to N, at time step N
  for j = (-N+1):N
    St[J[j]] = St[J[j]-1] * edx
  end

  ## initialise Array for the tree

  C = zeros(Float64, (2N+1, N+1))

  ## Option Values at expiry, time step = N, in rows -N to N

  for j = -N:N
    C[J[j], N+1] = max(0, St[J[j]] - K)
  end

  ## Stepping back through the lattice

  for i in range(N-1, -1, N)
    for j = -i:i
      C[J[j], i+1] = disc * (pu * C[J[j]+1, i+2] + pm * C[J[j], i+2] + pd * C[J[j]-1, i+2])
    end
  end

C
end

ext: jl


text:
"""
# Stipple

Stipple is a reactive UI library for Julia. It provides a rich API for building rich web UIs with 2-way bindings between
HTML UI elements and Julia. It requires minimum configuration, automatically setting up the WebSockets/Ajax communication
channels and automatically keeping the data in sync.

Stipple allows creating powerful reactive web data dashboards using only Julia coding. It employs a declarative
programming model, the framework taking care of the full data sync workflow.

Stipple uses Vue.js on the frontend and focuses mostly on the low level communication API, allowing data binding to HTML
elements, forms, etc. The core functionality is meant to be extended with specialized plugins, for example for powerful
UIs or plotting features. The plugins are easy to build, requiring a minimum Julia implementation that usually ties into
existing Vue.js libraries.
"""
module Stipple

using Logging, Reexport, Requires

@reexport using Observables
@reexport using Genie
@reexport using Genie.Renderer.Html
@reexport using JSON

include("NamedTuples.jl")
using .NamedTuples

const JSONParser = JSON
export JSONParser

# support for handling JS `undefined` values
export Undefined, UNDEFINED

struct Undefined
end

const UNDEFINED = Undefined()
const UNDEFINED_PLACEHOLDER = "__undefined__"
const UNDEFINED_VALUE = "undefined"

JSON.lower(x::Undefined) = UNDEFINED_PLACEHOLDER
Base.show(io::IO, x::Undefined) = Base.print(io, UNDEFINED_VALUE)

function Genie.Renderer.Html.attrparser(k::Symbol, v::JSONParser.JSONText) :: String
  if startswith(v.s, ":")
    ":$(k |> Genie.Renderer.Html.parseattr)=$(v.s[2:end]) "
  else
    "$(k |> Genie.Renderer.Html.parseattr)=$(v.s) "
  end
end

mutable struct Reactive{T} <: Observables.AbstractObservable{T}
  o::Observables.Observable{T}
  r_mode::Int
  no_backend_watcher::Bool
  no_frontend_watcher::Bool
  Reactive{T}() where {T} = new{T}(Observable{T}(), 0, false, false)
  Reactive{T}(o, no_bw::Bool = false, no_fw::Bool = false) where {T} = new{T}(o, 0, no_bw, no_fw)
  Reactive{T}(o, mode::Int, no_bw::Bool = false, no_fw::Bool = false) where {T} = new{T}(o, mode, no_bw, no_fw)
  Reactive{T}(o, mode::Int, updatemode::Int) where {T} = new{T}(o, mode, updatemode & NO_BACKEND_WATCHER != 0, updatemode & NO_FRONTEND_WATCHER != 0)

  # Construct an Reactive{Any} without runtime dispatch
  Reactive{Any}(@nospecialize(o)) = new{Any}(Observable{Any}(o), 0, false, false)
end

Reactive(r::T, arg1, args...) where T = convert(Reactive{T}, (r, arg1, args...))
Reactive(r::T) where T = convert(Reactive{T}, r)

Base.convert(::Type{T}, x::T) where {T<:Reactive} = x  # resolves ambiguity with convert(::Type{T}, x::T) in base/essentials.jl
Base.convert(::Type{T}, x) where {T<:Reactive} = T(x)

Base.convert(::Type{Reactive{T}}, (r, m)::Tuple{T, Int}) where T = m < 16 ? Reactive{T}(Observable(r), m, 0) : Reactive{T}(Observable(r), 0, m)
Base.convert(::Type{Reactive{T}}, (r, w)::Tuple{T, Bool}) where T = Reactive{T}(Observable(r), 0, w, false)
Base.convert(::Type{Reactive{T}}, (r, m, nw)::Tuple{T, Int, Bool}) where T = Reactive{T}(Observable(r), m, nw, false)
Base.convert(::Type{Reactive{T}}, (r, nbw, nfw)::Tuple{T, Bool, Bool}) where T = Reactive{T}(Observable(r), 0, nbw, nfw)
Base.convert(::Type{Reactive{T}}, (r, m, nbw, nfw)::Tuple{T, Int, Bool, Bool}) where T = Reactive{T}(Observable(r), m, nbw, nfw)
Base.convert(::Type{Reactive{T}}, (r, m, u)::Tuple{T, Int, Int}) where T = Reactive{T}(Observable(r), m, u)
Base.convert(::Type{Observable{T}}, r::Reactive{T}) where T = getfield(r, :o)

Base.getindex(r::Reactive{T}) where T = Base.getindex(getfield(r, :o))
Base.setindex!(r::Reactive{T}) where T = Base.setindex!(getfield(r, :o))

# pass indexing and property methods to referenced variable
function Base.getindex(r::Reactive{T}, arg1, args...) where T
  getindex(getfield(r, :o).val, arg1, args...)
end

function Base.setindex!(r::Reactive{T}, val, arg1, args...) where T
  setindex!(getfield(r, :o).val, val, arg1, args...)
  Observables.notify!(r)
end

Base.setindex!(r::Reactive, val, ::typeof(!)) = getfield(r, :o).val = val
Base.getindex(r::Reactive, ::typeof(!)) = getfield(r, :o).val

function Base.getproperty(r::Reactive{T}, field::Symbol) where T
  if field in (:o, :r_mode, :no_backend_watcher, :no_frontend_watcher) # fieldnames(Reactive)
    getfield(r, field)
  else
    if field == :val
      @warn """Reactive API has changed, use "[]" instead of ".val"!"""
      getfield(r, :o).val
    else
      getproperty(getfield(r, :o).val, field)
    end
  end
end

function Base.setproperty!(r::Reactive{T}, field::Symbol, val) where T
  if field in fieldnames(Reactive)
    setfield!(r, field, val)
  else
    if field == :val
      @warn """Reactive API has changed, use "setfield_withoutwatchers!() or o.val" instead of ".val"!"""
      getfield(r, :o).val = val
    else
      setproperty!(getfield(r, :o).val, field, val)
      Observables.notify!(r)
    end
  end
end

Observables.observe(r::Reactive{T}, args...; kwargs...) where T = Observables.observe(getfield(r, :o), args...; kwargs...)
Observables.listeners(r::Reactive{T}, args...; kwargs...) where T = Observables.listeners(getfield(r, :o), args...; kwargs...)

@static if isdefined(Observables, :appendinputs!)
    Observables.appendinputs!(r::Reactive{T}, obsfuncs) where T = Observables.appendinputs!(getfield(r, :o), obsfuncs)
end

# workaround for `map!()`, as long as Observables.MapUpdater is not patched to handle AbstractObservables (probably Observables <= v0.4.0)
import Base.map!
@inline Base.map!(f::F, r::Reactive, os...; update::Bool=true) where F = Base.map!(f::F, getfield(r, :o), os...; update=update)

const R = Reactive
const PUBLIC = 1
const PRIVATE = 2
const READONLY = 4
const JSFUNCTION = 8
const NO_BACKEND_WATCHER = 16
const NO_FRONTEND_WATCHER = 32
const NO_WATCHER = 48

OptDict = Dict{Symbol, Any}
opts(;kwargs...) = OptDict(kwargs...)

WEB_TRANSPORT = Genie.WebChannels

export R, Reactive, ReactiveModel, @R_str, @js_str, client_data
export PRIVATE, PUBLIC, READONLY, JSFUNCTION, NO_WATCHER, NO_BACKEND_WATCHER, NO_FRONTEND_WATCHER
export newapp
export onbutton
export @kwredef

#===#

function __init__()
  Genie.config.websockets_server = true
  @require OffsetArrays  = "6fe1bfb0-de20-5000-8ca7-80f57d26f881" function convertvalue(targetfield::Union{Ref{T}, Reactive{T}}, value) where T <: OffsetArrays.OffsetArray
    valtype = eltype(targetfield)
    a = if value isa AbstractArray
      convert(Array{eltype(valtype)}, value)
    else
      eltype(valtype)[value]
    end
    if ! isa(value, OffsetArrays.OffsetArray)
      o = targetfield[].offsets
      OffsetArrays.OffsetArray(a, OffsetArrays.Origin(1 .+ o))
    else
      value
    end
  end
end

#===#

"""
    `type ReactiveModel`

The abstract type that is inherited by Stipple models. Stipple models are used for automatic 2-way data sync and data
exchange between the Julia backend and the JavaScript/Vue.js frontend.

### Example

```julia
Base.@kwdef mutable struct HelloPie <: ReactiveModel
  plot_options::R{PlotOptions} = PlotOptions(chart_type=:pie, chart_width=380, chart_animations_enabled=true,
                                            stroke_show = false, labels=["Slice A", "Slice B"])
  piechart::R{Vector{Int}} = [44, 55]
  values::R{String} = join(piechart, ",")
end
```
"""
abstract type ReactiveModel end
mutable struct Settings
  readonly_pattern
  private_pattern
end
Settings(; readonly_pattern = r"_$", private_pattern = r"__$") = Settings(readonly_pattern, private_pattern)

#===#

"""
    `const JS_SCRIPT_NAME`

The name of the dynamically generated JavaScript file used for data sync.
"""
const JS_SCRIPT_NAME = "stipple.js"

"""
    `const JS_DEBOUNCE_TIME`

Debounce time used to indicate the minimum frequency for sending data payloads to the backend (for example to batch send
payloads when the user types into an text field, to avoid overloading the server).
"""
const JS_DEBOUNCE_TIME = 300 #ms
const SETTINGS = Settings()

#===#

"""
    `function render`

Abstract function. Needs to be specialized by plugins. It is automatically invoked by `Stipple` to serialize a Julia
data type (corresponding to the fields in the `ReactiveModel` instance) to JavaScript/JSON. In general the specialized
methods should return a Julia `Dict` which are automatically JSON encoded by `Stipple`. If custom JSON serialization is
required for certain types in the resulting `Dict`, specialize `JSON.lower` for that specific type.

### Example

```julia
function Stipple.render(ps::PlotSeries, fieldname::Union{Symbol,Nothing} = nothing)
  Dict(:name => ps.name, ps.plotdata.key => ps.plotdata.data)
end
```

#### Specialized JSON rendering for `Undefined`

```julia
JSON.lower(x::Undefined) = "__undefined__"
```
"""
function render end

"""
    `function update! :: {M<:ReactiveModel}`

Abstract function used to update the values of the fields in the `ReactiveModel` based on the data from the frontend.
Can be specialized for dedicated types, but it is usually not necessary. If specialized, it must return the update
instance of `ReactiveModel` provided as the first parameter.

### Example

```julia
function update!(model::M, field::Any, newval::T, oldval::T)::M where {T,M<:ReactiveModel}
  setfield!(model, field, newval)

  model
end
````
"""
function update! end

"""
    `function watch`

Abstract function. Can be used by plugins to define custom Vue.js watch functions.
"""
function watch end

"""
    `function js_methods(app::T) where {T<:ReactiveModel}`

Defines js functions for the `methods` section of the vue element.

### Example

```julia
js_methods(app::MyDashboard) = \"\"\"
  mysquare: function (x) {
    return x^2
  }
  myadd: function (x, y) {
    return x + y
  }
\"\"\"
```
"""
function js_methods(app::T)::String where {T<:ReactiveModel}
  ""
end

"""
    `function js_computed(app::T) where {T<:ReactiveModel}`

Defines js functions for the `computed` section of the vue element.
These properties are updated every time on of the inner parameters changes its value.

### Example

```julia
js_computed(app::MyDashboard) = \"\"\"
  fullName: function () {
    return this.firstName + ' ' + this.lastName
  }
\"\"\"
```
"""
function js_computed(app::T)::String where {T<:ReactiveModel}
  ""
end

const jscomputed = js_computed

"""
    `function js_watch(app::T) where {T<:ReactiveModel}`

Defines js functions for the `watch` section of the vue element.
These functions are called every time the respective property changes.

### Example

Updates the `fullName` every time `firstName` or `lastName` changes.

```julia
js_watch(app::MyDashboard) = \"\"\"
  firstName: function (val) {
    this.fullName = val + ' ' + this.lastName
  },
  lastName: function (val) {
    this.fullName = this.firstName + ' ' + val
  }
\"\"\"
```
"""
function js_watch(m::T)::String where {T<:ReactiveModel}
  ""
end

const jswatch = js_watch

"""
    `function js_created(app::T)::String where {T<:ReactiveModel}`

Defines js statements for the `created` section of the vue element.
They are executed directly after the creation of the vue element.

### Example

```julia
js_created(app::MyDashboard) = \"\"\"
    if (this.cameraon) { startcamera() }
\"\"\"
```
"""
function js_created(app::T)::String where {T<:ReactiveModel}
  ""
end

const jscreated = js_created

"""
    `function js_mounted(app::T)::String where {T<:ReactiveModel}`

Defines js statements for the `mounted` section of the vue element.
They are executed directly after the mounting of the vue element.

### Example

```julia
js_created(app::MyDashboard) = \"\"\"
    if (this.cameraon) { startcamera() }
\"\"\"
```
"""
function js_mounted(app::T)::String where {T<:ReactiveModel}
  ""
end

const jsmounted = js_mounted

"""
    `function client_data(app::T)::String where {T<:ReactiveModel}`

Defines additional data that will only be visible by the browser.

It is meant to keep volatile data, e.g. form data that needs to pass a validation first.
In order to use the data you most probably also want to define [`js_methods`](@ref)
### Example

```julia
import Stipple.client_data
client_data(m::Example) = client_data(client_name = js"null", client_age = js"null", accept = false)
```
will define the additional fields `client_name`, `clientage` and `accept` for the model `Example`. These should, of course, not overlap with existing fields of your model.
"""
client_data(app::T) where T <: ReactiveModel = Dict{String, Any}()

client_data(;kwargs...) = Dict{String, Any}([String(k) => v for (k, v) in kwargs]...)

#===#

const COMPONENTS = Dict()

"""
    `function register_components(model::Type{M}, keysvals::AbstractVector) where {M<:ReactiveModel}`

Utility function for adding Vue components that need to be registered with the Vue.js app.
This is usually needed for registering components provided by Stipple plugins.

### Example

```julia
Stipple.register_components(HelloPie, StippleCharts.COMPONENTS)
```
"""
function register_components(model::Type{M}, keysvals::AbstractVector) where {M<:ReactiveModel}
  haskey(COMPONENTS, model) || (COMPONENTS[model] = Any[])
  push!(COMPONENTS[model], keysvals...)
end

"""
    `function components(m::Type{M})::String where {M<:ReactiveModel}`
    `function components(app::M)::String where {M<:ReactiveModel}`

JSON representation of the Vue.js components registered for the `ReactiveModel` `M`.
"""
function components(m::Type{M})::String where {M<:ReactiveModel}
  haskey(COMPONENTS, m) || return ""

  replace(Dict(COMPONENTS[m]...) |> JSON.json, "\""=>"") |> string
end

function components(app::M)::String where {M<:ReactiveModel}
  components(M)
end

#===#

"""
    `setindex_withoutwatchers!(field::Reactive, val; notify=(x)->true)`
    `setindex_withoutwatchers!(field::Reactive, val, keys::Int...; notify=(x)->true)`

Change the content of a Reactive field without triggering the listeners.
If keys are specified, only these listeners are exempted from triggering.
"""
function setindex_withoutwatchers!(field::Reactive{T}, val, keys::Int...; notify=(x)->true) where T
  count = 1
  field.o.val = val
  length(keys) == 0 && return field

  for f in Observables.listeners(field.o)
    if in(count, keys)
      count += 1

      continue
    end

    if notify(f)
      try
        f(val)
      catch ex
        isa(ex, MethodError) && string(ex.f) == string(f) ? Base.invokelatest(f, val) : rethrow(ex)
      end
    end

    count += 1
  end

  return field
end

function Base.setindex!(r::Reactive{T}, val, args::Vector{Int}; notify=(x)->true) where T
  setindex_withoutwatchers!(r, val, args...)
end

"""
    `setfield_withoutwatchers!(app::ReactiveModel, field::Symmbol, val; notify=(x)->true)``
    `setfield_withoutwatchers!(app::ReactiveModel, field::Symmbol, val, keys...; notify=(x)->true)`

Change the field of a ReactiveModel without triggering the listeners.
If keys are specified, only these listeners are exempted from triggering.
"""
function setfield_withoutwatchers!(app::T, field::Symbol, val, keys...; notify=(x)->true) where T <: ReactiveModel
  f = getfield(app, field)
  if f isa Reactive
    setindex_withoutwatchers!(f, val, keys...; notify = notify)
  else
    setfield!(app, field, val)
  end
  return app
end

#===#

include("Typography.jl")
include("Elements.jl")
include("Layout.jl")

@reexport using .Typography
@reexport using .Elements
@reexport using .Layout

#===#

function convertvalue(targetfield::Any, value)
  valtype = eltype(targetfield)
  try
    if valtype <: AbstractFloat && typeof(value) <: Integer
      convert(valtype, value)
    elseif valtype <: AbstractArray
      a = if value isa AbstractArray
        convert(Array{eltype(valtype)}, value)
      else
        eltype(valtype)[value]
      end
    else
      Base.parse(valtype, value)
    end
  catch ex
    @error ex
    value
  end
end

"""
    `function update!(model::M, field::Symbol, newval::T, oldval::T)::M where {T,M<:ReactiveModel}`
    `function update!(model::M, field::Reactive, newval::T, oldval::T)::M where {T,M<:ReactiveModel}`
    `function update!(model::M, field::Any, newval::T, oldval::T)::M where {T,M<:ReactiveModel}`

Sets the value of `model.field` from `oldval` to `newval`. Returns the upated `model` instance.
"""
function update!(model::M, field::Symbol, newval::T1, oldval::T2)::M where {T1, T2, M<:ReactiveModel}
  f = getfield(model, field)
  if f isa Reactive
    f.r_mode == PRIVATE || f.no_backend_watcher ? f[] = newval : setindex_withoutwatchers!(f, newval, 1)
  else
    setfield!(model, field, newval)
  end
  model
end

function update!(model::M, field::Reactive{T}, newval::T, oldval::T)::M where {T, M<:ReactiveModel}
  field.r_mode == PRIVATE || field.no_backend_watcher ? field[] = newval : setindex_withoutwatchers!(field, newval, 1)

  model
end

function update!(model::M, field::Any, newval::T, oldval::T)::M where {T,M<:ReactiveModel}
  setfield!(model, field, newval)

  model
end

#===#

"""
    `function watch(vue_app_name::String, fieldtype::Any, fieldname::Symbol, channel::String, debounce::Int, model::M)::String where {M<:ReactiveModel}`

Sets up default Vue.js watchers so that when the value `fieldname` of type `fieldtype` in model `vue_app_name` is
changed on the frontend, it is pushed over to the backend using `channel`, at a `debounce` minimum time interval.
"""
function watch(vue_app_name::String, fieldname::Symbol, channel::String, debounce::Int, model::M)::String where {M<:ReactiveModel}
  js_channel = channel == "" ? "window.Genie.Settings.webchannels_default_route" : "'$channel'"

  output = """
    $vue_app_name.\$watch(function(){return this.$fieldname}, _.debounce(function(newVal, oldVal){
      Genie.WebChannels.sendMessageTo($js_channel, 'watchers', {'payload': {'field':'$fieldname', 'newval': newVal, 'oldval': oldVal}});
    }, $debounce), {deep: true});
  """
  # in production mode vue does not fill `this.expression` in the watcher, so we do it manually
  if Genie.Configuration.isprod()
    output *= "$vue_app_name._watchers[$vue_app_name._watchers.length - 1].expression = 'function(){return this.$fieldname}'"
  end

  output *= "\n\n"
end

#===#

function Base.parse(::Type{T}, v::T) where {T}
  v::T
end

"""
    `function init(model::M, ui::Union{String,Vector} = ""; vue_app_name::String = Stipple.Elements.root(model),
                    endpoint::String = JS_SCRIPT_NAME, channel::String = Genie.config.webchannels_default_route,
                    debounce::Int = JS_DEBOUNCE_TIME, transport::Module = Genie.WebChannels)::M where {M<:ReactiveModel}`

Initializes the reactivity of the model `M` by setting up the custom JavaScript for integrating with the Vue.js
frontend and perform the 2-way backend-frontend data sync. Returns the instance of the model.

### Example

```julia
hs_model = Stipple.init(HelloPie())
```
"""
function init(model::M, ui::Union{String,Vector} = ""; vue_app_name::String = Stipple.Elements.root(model),
              endpoint::String = JS_SCRIPT_NAME, channel::String = Genie.config.webchannels_default_route,
              debounce::Int = JS_DEBOUNCE_TIME, transport::Module = Genie.WebChannels,
              parse_errors::Bool = false)::M where {M<:ReactiveModel}

  global WEB_TRANSPORT = transport
  transport == Genie.WebChannels || (Genie.config.websockets_server = false)

  deps_routes(channel)

  Genie.Router.channel("/$(channel)/watchers") do
    payload = Genie.Requests.payload(:payload)["payload"]
    client = Genie.Requests.wsclient()

    # if only elements of the array change, oldval and newval are identical
    ! isa(payload["newval"], Array) && ! isa(payload["newval"], Dict) && payload["newval"] == payload["oldval"] && return "OK"

    field = Symbol(payload["field"])

    #check if field exists
    hasfield(M, field) || return "OK"
    valtype = Dict(zip(fieldnames(M), M.types))[field]
    val = valtype <: Reactive ? getfield(model, field) : Ref{valtype}(getfield(model, field))

    # reject non-public types
    if val isa Reactive
      val.r_mode == PUBLIC || return "OK"
    elseif occursin(SETTINGS.readonly_pattern, String(field)) || occursin(SETTINGS.private_pattern, String(field))
      return "OK"
    end

    newval = convertvalue(val, payload["newval"])
    oldval = convertvalue(val, payload["oldval"])

    push!(model, field => newval, channel = channel, except = client)
    update!(model, field, newval, oldval)

    "OK"
  end

  ep = channel == Genie.config.webchannels_default_route ? endpoint : "js/$channel/$endpoint"
  Genie.Router.route("/$(ep)") do
    Stipple.Elements.vue_integration(model, vue_app_name = vue_app_name, endpoint = ep, channel = "", debounce = debounce) |> Genie.Renderer.Js.js
  end

  setup(model, channel)
end

"""
    `function setup(model::M, channel = Genie.config.webchannels_default_route)::M where {M<:ReactiveModel}`

Configures the reactive handlers for the reactive properties of the model. Called internally.
"""
function setup(model::M, channel = Genie.config.webchannels_default_route)::M where {M<:ReactiveModel}
  for field in fieldnames(M)
    f = getproperty(model, field)

    isa(f, Reactive) || continue

    if f.r_mode == 0
      if occursin(SETTINGS.private_pattern, String(field))
        f.r_mode = PRIVATE
      elseif occursin(SETTINGS.readonly_pattern, String(field))
        f.r_mode = READONLY
      else
        f.r_mode = PUBLIC
      end
    end
    f.r_mode == PRIVATE || f.no_backend_watcher && continue

    on(f) do _
      push!(model, field => f, channel = channel)
    end
  end

  model
end

#===#

"""
    `Base.push!(app::M, vals::Pair{Symbol,T}; channel::String,
                except::Union{Genie.WebChannels.HTTP.WebSockets.WebSocket,Nothing,UInt}) where {T,M<:ReactiveModel}`

Pushes data payloads over to the frontend by broadcasting the `vals` through the `channel`.
"""
function Base.push!(app::M, vals::Pair{Symbol,T};
                    channel::String = Genie.config.webchannels_default_route,
                    except::Union{Genie.WebChannels.HTTP.WebSockets.WebSocket,Nothing,UInt} = nothing) where {T,M<:ReactiveModel}
  try
    WEB_TRANSPORT.broadcast(channel,
                                    JSON.json(Dict( "key" => julia_to_vue(vals[1]),
                                                    "value" => Stipple.render(vals[2], vals[1]))),
                                    except = except)
  catch
  end
end

function Base.push!(app::M, vals::Pair{Symbol,Reactive{T}};
                    channel::String = Genie.config.webchannels_default_route,
                    except::Union{Genie.WebChannels.HTTP.WebSockets.WebSocket,Nothing,UInt} = nothing) where {T,M<:ReactiveModel}
                    v = vals[2].r_mode != JSFUNCTION ? vals[2][] : replace_jsfunction(vals[2][])
  push!(app, Symbol(julia_to_vue(vals[1])) => v, channel = channel, except = except)
end

#===#

RENDERING_MAPPINGS = Dict{String,String}()
mapping_keys() = collect(keys(RENDERING_MAPPINGS))

"""
    `function rendering_mappings(mappings = Dict{String,String})`

Registers additional `mappings` as Julia to Vue properties mappings  (eg `foobar` to `foo-bar`).
"""
function rendering_mappings(mappings = Dict{String,String})
  merge!(RENDERING_MAPPINGS, mappings)
end

"""
    `function julia_to_vue(field, mapping_keys = mapping_keys())`

Converts Julia names to Vue names (eg `foobar` to `foo-bar`).
"""
function julia_to_vue(field, mapping_keys = mapping_keys()) :: String
  if in(string(field), mapping_keys)
    parts = split(RENDERING_MAPPINGS[string(field)], "-")

    if length(parts) > 1
      extraparts = map((x) -> uppercasefirst(string(x)), parts[2:end])
      string(parts[1], join(extraparts))
    else
      parts |> string
    end
  else
    field |> string
  end
end

"""
    `function Stipple.render(app::M, fieldname::Union{Symbol,Nothing} = nothing)::Dict{Symbol,Any} where {M<:ReactiveModel}`

Renders the Julia `ReactiveModel` `app` as the corresponding Vue.js JavaScript code.
"""
function Stipple.render(app::M, fieldname::Union{Symbol,Nothing} = nothing)::Dict{Symbol,Any} where {M<:ReactiveModel}
  result = Dict{String,Any}()
  for field in fieldnames(typeof(app))
    f = getfield(app, field)

    occursin(SETTINGS.private_pattern, String(field)) && continue
    f isa Reactive && f.r_mode == PRIVATE && continue

    result[julia_to_vue(field)] = Stipple.render(f, field)
  end

  vue = Dict(:el => Elements.elem(app), :mixins => JSONText("[watcherMixin, reviveMixin]"), :data => merge(result, client_data(app)))

  isempty(components(app)   |> strip)   || push!(vue, :components => components(app))
  isempty(js_methods(app)   |> strip)   || push!(vue, :methods    => JSONText("{ $(js_methods(app)) }"))
  isempty(js_computed(app)  |> strip)   || push!(vue, :computed   => JSONText("{ $(js_computed(app)) }"))
  isempty(js_watch(app)     |> strip)   || push!(vue, :watch      => JSONText("{ $(js_watch(app)) }"))
  isempty(js_created(app)   |> strip)   || push!(vue, :created    => JSONText("function(){ $(js_created(app)); }"))
  isempty(js_mounted(app)   |> strip)   || push!(vue, :mounted    => JSONText("function(){ $(js_mounted(app)); }"))

  vue
end

"""
    `function Stipple.render(val::T, fieldname::Union{Symbol,Nothing} = nothing) where {T}`

Default rendering of value types. Specialize `Stipple.render` to define custom rendering for your types.
"""
function Stipple.render(val::T, fieldname::Union{Symbol,Nothing} = nothing) where {T}
  val
end

"""
    `function Stipple.render(o::Reactive{T}, fieldname::Union{Symbol,Nothing} = nothing) where {T}`

Default rendering of `Reactive` values. Specialize `Stipple.render` to define custom rendering for your type `T`.
"""
function Stipple.render(o::Reactive{T}, fieldname::Union{Symbol,Nothing} = nothing) where {T}
  Stipple.render(o[], fieldname)
end

"""
```
function parse_jsfunction(s::AbstractString)
```
Checks whether the string is a valid js function and returns a `Dict` from which a reviver function
in the backend can construct a function.
"""
function parse_jsfunction(s::AbstractString)
    # look for classical function definition
    m = match( r"function\s*\(([^)]*)\)\s*{(.*)}", s)
    !isnothing(m) && length(m.captures) == 2 && return opts(arguments=m[1], body=m[2])

    # look for pure function definition
    m = match( r"\s*\(?([^=)]*?)\)?\s*=>\s*({*.*?}*)\s*$" , s )
    (isnothing(m) || length(m.captures) != 2) && return nothing

    # if pure function body is without curly brackets, add a `return`, otherwise strip the brackets
    # Note: for utf-8 strings m[2][2:end-1] will fail if the string ends with a wide character, e.g. ϕ
    body = startswith(m[2], "{") ? m[2][2:prevind(m[2], lastindex(m[2]))] : "return " * m[2]
    return opts(arguments=m[1], body=body)
end

"""
```
function replace_jsfunction!(js::Union{Dict, JSONText})
```
Replaces all JSONText values that contain a valid js function by a `Dict` that codes the function for a reviver.
For JSONText variables it encapsulates the dict in a JSONText to make the function type stable.
"""
# fallback is identity function
replace_jsfunction!(x) = x

function replace_jsfunction!(d::Dict)
    for (k,v) in d
        if isa(v, Dict) || isa(v, Array)
            replace_jsfunction!(v)
        elseif isa(v, JSONText)
            jsfunc = parse_jsfunction(v.s)
            isnothing(jsfunc) || ( d[k] = opts(jsfunction=jsfunc) )
        end
    end
    return d
end

function replace_jsfunction!(v::Array)
  replace_jsfunction!.(v)
end

function replace_jsfunction(d::Dict)
  replace_jsfunction!(deepcopy(d))
end

function replace_jsfunction(v::Vector)
  replace_jsfunction!.(deepcopy(v))
end

function replace_jsfunction(js::JSONText)
    jsfunc = parse_jsfunction(js.s)
    isnothing(jsfunc) ? js : JSONText(json(opts(jsfunction=jsfunc)))
end
#===#

const DEPS = Function[]

"""
    `function deps_routes(channel::String = Genie.config.webchannels_default_route) :: Nothing`

Registers the `routes` for all the required JavaScript dependencies (scripts).
"""
function deps_routes(channel::String = Genie.config.webchannels_default_route) :: Nothing
  VUEJS = Genie.Configuration.isprod() ? "vue.min.js" : "vue.js"

  Genie.Router.route("/js/stipple/$VUEJS") do
    Genie.Renderer.WebRenderable(
      read(joinpath(@__DIR__, "..", "files", "js", VUEJS), String),
      :javascript) |> Genie.Renderer.respond
  end

  Genie.Router.route("/js/stipple/vue_filters.js") do
    Genie.Renderer.WebRenderable(
      read(joinpath(@__DIR__, "..", "files", "js", "vue_filters.js"), String),
      :javascript) |> Genie.Renderer.respond
  end

  Genie.Router.route("/js/stipple/underscore-min.js") do
    Genie.Renderer.WebRenderable(
      read(joinpath(@__DIR__, "..", "files", "js", "underscore-min.js"), String),
      :javascript) |> Genie.Renderer.respond
  end

  Genie.Router.route("/js/stipple/stipplecore.js") do
    Genie.Renderer.WebRenderable(
      read(joinpath(@__DIR__, "..", "files", "js", "stipplecore.js"), String),
      :javascript) |> Genie.Renderer.respond
  end

  (WEB_TRANSPORT == Genie.WebChannels ? Genie.Assets.channels_support(channel) : Genie.Assets.webthreads_support(channel))

  nothing
end

"""
    `function deps(channel::String = Genie.config.webchannels_default_route) :: String`

Outputs the HTML code necessary for injecting the dependencies in the page (the <script> tags).
"""
function deps(channel::String = Genie.config.webchannels_default_route) :: String
  VUEJS = Genie.Configuration.isprod() ? "vue.min.js" : "vue.js"

  endpoint = (channel == Genie.config.webchannels_default_route) ?
              Stipple.JS_SCRIPT_NAME :
              "js/$(channel)/$(Stipple.JS_SCRIPT_NAME)"

  string(
    (WEB_TRANSPORT == Genie.WebChannels ? Genie.Assets.channels_support(channel) : Genie.Assets.webthreads_support(channel)),
    Genie.Renderer.Html.script(src="$(Genie.config.base_path)js/stipple/underscore-min.js"),
    Genie.Renderer.Html.script(src="$(Genie.config.base_path)js/stipple/$VUEJS"),
    join([f() for f in DEPS], "\n"),
    Genie.Renderer.Html.script(src="$(Genie.config.base_path)js/stipple/stipplecore.js", defer=true),
    Genie.Renderer.Html.script(src="$(Genie.config.base_path)js/stipple/vue_filters.js", defer=true),

    # if the model is not configured and we don't generate the stipple.js file, no point in requesting it
    in(Symbol("get_$(replace(endpoint, '/' => '_'))"), Genie.Router.named_routes() |> keys |> collect)
      ?
      Genie.Renderer.Html.script(src="$(Genie.config.base_path)$(endpoint)?v=$(Genie.Configuration.isdev() ? rand() : 1)",
                                    defer=true, onload="Stipple.init({theme: 'stipple-blue'});")
      :
      begin
        @warn "The Reactive Model is not initialized - make sure you call Stipple.init(YourModel()) to initialize it"
        ""
      end
  )
end

macro R_str(s)
  :(Symbol($s))
end

#== ==#

# add a method to Observables.on to accept inverted order of arguments similar to route()
import Observables.on
on(observable::Observables.AbstractObservable, f::Function; weak = false) = on(f, observable; weak = weak)

"""
    `onbutton(f::Function, button::R{Bool}; async = false, weak = false)`

Links a function to a reactive boolean parameter, typically a representing a button of an app.
After the function is called, the parameter is set back to false. The `async` keyword
specifies whether the call should be made asynchroneously or not.

### Example

```julia
onbutton(model.save_button) do
  # save what has to be saved
end
```
"""
onbutton(f::Function, button::R{Bool}; async = false, kwargs...) = on(button; kwargs...) do pressed
  pressed || return
  if async
      @async begin
          try
            f()
          catch ex
            warn(ex)
          end
          button[] = false
      end
  else
      try
        f()
      catch ex
        @warn(ex)
      end
      button[] = false
  end
  return
end

onbutton(button::R{Bool}, f::Function; kwargs...) = onbutton(f, button; kwargs...)

"""
    `@js_str -> JSONText`

Construct a JSONText, such as `js"button=false"`, without interpolation and unescaping
(except for quotation marks `"`` which still has to be escaped). Avoiding escaping `"`` can be done by
`js\"\"\"alert("Hello World")\"\"\"`.
"""
macro js_str(expr)
  :( JSONText($(esc(expr))) )
end

"""
    `@kwredef(expr)`

Helper function during development that is a one-to-one replacement for `@kwdef` but allows for redefinition of the struct.

Internally it defines a new struct with a number appended to the original struct name
and assigns this struct to a variable with the original struct name.
"""
macro kwredef(expr)
  expr = macroexpand(__module__, expr) # to expand @static
  expr isa Expr && expr.head === :struct || error("Invalid usage of @kwredef")
  expr = expr::Expr

  t = expr.args; n = 2
  if t[n] isa Expr && t[n].head === :<:
      t = t[n].args
      n = 1
  end
  curly = t[n] isa Expr && t[n].head === :curly
  if curly
      t = t[n].args
      n=1
  end

  T_old = t[n]
  t[n] = T_new = gensym(T_old)

  esc(quote
    Base.@kwdef $expr
    $T_old = $T_new
    $curly ? $T_new.body.name.name = $(QuoteNode(T_old)) : $T_new.name.name = $(QuoteNode(T_old)) # fix the name
  end)
end

"""
    `Stipple.@kwdef`

Helper function for model definition that acts as a one-to-one replacement for `Base.@kwdef`.

When `Genie.Configuration.isprod() == true` this macro calls `@kwredef` and allows for redefinition of models.
Otherwise it calls `Base.@kwdef`.
"""
macro kwdef(expr)
  esc(quote
    Genie.Configuration.isprod() ? Base.@kwdef($expr) : Stipple.@kwredef($expr)
  end)
end

end


text:
#= src/trapezoidal.jl
=#

abstract type Trapezoidal{T} <: AbstractMatrix{T} end

mutable struct LowerTrapezoidal{T} <: Trapezoidal{T}
  data::Matrix{T}
end

mutable struct UpperTrapezoidal{T} <: Trapezoidal{T}
  data::Matrix{T}
end

LowerTrapezoidal(data::AbstractMatrix) = LowerTrapezoidal(Matrix(data))
UpperTrapezoidal(data::AbstractMatrix) = UpperTrapezoidal(Matrix(data))

conj!(A::LowerTrapezoidal) = LowerTrapezoidal(conj!(A.data))
conj!(A::UpperTrapezoidal) = UpperTrapezoidal(conj!(A.data))
conj(A::Trapezoidal) = conj!(copy(A))

convert(::Type{LowerTrapezoidal{T}}, A::LowerTrapezoidal) where {T} =
  LowerTrapezoidal(convert(Array{T}, A.data))
convert(::Type{UpperTrapezoidal{T}}, A::UpperTrapezoidal) where {T} =
  UpperTrapezoidal(convert(Array{T}, A.data))
convert(::Type{Trapezoidal{T}}, A::LowerTrapezoidal) where {T} =
  convert(LowerTrapezoidal{T}, A)
convert(::Type{Trapezoidal{T}}, A::UpperTrapezoidal) where {T} =
  convert(UpperTrapezoidal{T}, A)
convert(::Type{Array}, A::Trapezoidal) = Matrix(A)
convert(::Type{Array{T}}, A::Trapezoidal) where {T} = convert(Array{T}, Matrix(A))
convert(::Type{Matrix}, A::Trapezoidal) = Matrix(A)
convert(::Type{Matrix{T}}, A::Trapezoidal) where {T} = convert(Array{T}, Matrix(A))

copy(A::LowerTrapezoidal) = LowerTrapezoidal(copy(A.data))
copy(A::UpperTrapezoidal) = UpperTrapezoidal(copy(A.data))

adjoint(A::LowerTrapezoidal) = UpperTrapezoidal(A.data')
adjoint(A::UpperTrapezoidal) = LowerTrapezoidal(A.data')
transpose(A::LowerTrapezoidal) = UpperTrapezoidal(transpose(A.data))
transpose(A::UpperTrapezoidal) = LowerTrapezoidal(transpose(A.data))

Matrix(A::LowerTrapezoidal) = tril(A.data)
Matrix(A::UpperTrapezoidal) = triu(A.data)

getindex(A::LowerTrapezoidal{T}, i::Integer, j::Integer) where {T} =
  i >= j ? A.data[i,j] : zero(T)
getindex(A::UpperTrapezoidal{T}, i::Integer, j::Integer) where {T} =
  i <= j ? A.data[i,j] : zero(T)

imag(A::LowerTrapezoidal) = LowerTrapezoidal(imag(A.data))
imag(A::UpperTrapezoidal) = UpperTrapezoidal(imag(A.data))

real(A::LowerTrapezoidal) = LowerTrapezoidal(real(A.data))
real(A::UpperTrapezoidal) = UpperTrapezoidal(real(A.data))

setindex!(A::LowerTrapezoidal, x, i::Integer, j::Integer) =
  i >= j ? (A.data[i,j] = x; A) : throw(BoundsError)
setindex!(A::UpperTrapezoidal, x, i::Integer, j::Integer) =
  i <= j ? (A.data[i,j] = x; A) : throw(BoundsError)

size(A::Trapezoidal, args...) = size(A.data, args...)

# BLAS/LAPACK routines

## LowerTrapezoidal left-multiplication
  function mul!(y::AbstractVector{T}, A::LowerTrapezoidal{T}, x::AbstractVector{T}) where T
    m, n = size(A)
    y[1:n] = x
    lmul!(LowerTriangular(view(A.data,1:n,:)), view(y,1:n))
    mul!(view(y,n+1:m), view(A.data,n+1:m,:), x)
    y
  end
  function mul!(C::AbstractMatrix{T}, A::LowerTrapezoidal{T}, B::AbstractMatrix{T}) where T
    m, n = size(A)
    C[1:n,:] = B
    lmul!(LowerTriangular(view(A.data,1:n,:)), view(C,1:n,:))
    mul!(view(C,n+1:m,:), view(A.data,n+1:m,:), B)
    C
  end
  for Adj in (:Adjoint, :Transpose)
      @eval function mul!(C::AbstractMatrix{T}, A::LowerTrapezoidal{T}, B::$Adj{T,<:AbstractMatrix{T}}) where T
        m, n = size(A)
        C[1:n,:] = B
        lmul!(LowerTriangular(view(A.data,1:n,:)), view(C,1:n,:))
        mul!(view(C,n+1:m,:), view(A.data,n+1:m,:), B)
        C
      end
  end
  for (Adj, g) in ((:Adjoint, :adjoint!), (:Transpose, :transpose!))
    @eval function mul!(C::AbstractMatrix{T}, A::$Adj{T,LowerTrapezoidal{T}}, B::AbstractMatrix{T}) where T
        m, n = size(A)
        $g(view(C,1:n,:), B)
        lmul!(LowerTriangular(view(A.data,1:n,:)), view(C,1:n,:))
        mul!(view(C,n+1:m,:), $Adj(view(A.data,n+1:m,:)), B)
        C
    end
  end

  for (Adj, trans) in ((:Adjoint, 'C'), (:Transpose, 'T'))
    @eval begin
      function mul!(y::AbstractVector{T}, adjA::$Adj{T,LowerTrapezoidal{T}}, x::AbstractVector{T}) where T<:BlasFloat
        A = parent(adjA)
        m, n = size(A)
        copyto!(y, view(x,1:n))
        BLAS.trmv!('L', $trans, 'N', view(A.data,1:n,:), y)
        BLAS.gemv!($trans, one(T), view(A.data,n+1:m,:), view(x,n+1:m), one(T), y)
      end
      function mul!(C::AbstractMatrix{T}, adjA::$Adj{T,LowerTrapezoidal{T}}, B::AbstractMatrix{T}) where T<:BlasFloat
        A = parent(adjA)
        m, n = size(A)
        copyto!(C, view(B,1:n,:))
        BLAS.trmm!('L', 'L', $trans, 'N', one(T), view(A.data,1:n,:), C)
        BLAS.gemm!($trans, 'N', one(T), view(A.data,n+1:m,:), view(B,n+1:m,:), one(T), C)
      end
    end
  end

  for (Adj, g, trans) in ((:Adjoint, :adjoint!, 'C'),
                        (:Transpose, :transpose!,  'T'))
    @eval begin
      function mul!(C::AbstractMatrix{T}, adjA::$Adj{T,LowerTrapezoidal{T}}, adjB::$Adj{T,<:AbstractMatrix{T}}) where T<:BlasFloat
        A = parent(adjA)
        B = parent(adjB)
        m, n = size(A)
        $g(C, view(B,:,1:n))
        BLAS.trmm!(
          'L', 'L', $trans, 'N', one(T), view(A.data,1:n,:), view(C,1:n,:))
        BLAS.gemm!(
          $trans, $trans, one(T), view(A.data,n+1:m,:), view(B,:,n+1:m), one(T), C)
      end
    end
  end

  ## LowerTrapezoidal right-multiplication

  function mul!(C::AbstractMatrix{T}, A::AbstractMatrix{T}, B::LowerTrapezoidal{T}) where T<:BlasFloat
    m, n = size(B)
    copyto!(C, view(A,:,1:n))
    rmul!(C, LowerTriangular(view(B.data,1:n,:)))
    BLAS.gemm!('N', 'N', one(T), view(A,:,n+1:m), view(B.data,n+1:m,:), one(T), C)
  end

  for (Adj, trans) in ((:Adjoint, 'C'), (:Transpose, 'T'))
    @eval begin
      function mul!(C::AbstractMatrix{T}, A::AbstractMatrix{T}, adjB::$Adj{T,LowerTrapezoidal{T}}) where T<:BlasFloat
        B = parent(adjB)
        m, n = size(B)
        copyto!(view(C,:,1:n), A)
        BLAS.trmm!(
          'R', 'L', $trans, 'N', one(T), view(B.data,1:n,:), view(C,:,1:n))
        mul!(view(C,:,n+1:m), A, $Adj(view(B.data,n+1:m,:)))
        C
      end
    end
  end

  for (Adj, g, trans) in ((:Adjoint, :adjoint!, 'C'),
                        (:Transpose, :transpose!,  'T'))
    @eval begin
      function mul!(
          C::AbstractMatrix{T}, adjA::$Adj{T,<:AbstractMatrix{T}}, B::LowerTrapezoidal{T}) where T<:BlasFloat
        A = parent(adjA)
        m, n = size(B)
        $g(C, view(A,1:n,:))
        BLAS.trmm!('R', 'L', 'N', 'N', one(T), view(B.data,1:n,:), C)
        BLAS.gemm!(
          $trans, 'N', one(T), view(A,n+1:m,:), view(B.data,n+1:m,:), one(T), C)
      end
    end
  end

  for (Adj, g, trans) in ((:Adjoint, :adjoint!, 'C'),
                        (:Transpose, :transpose!,  'T'))
    @eval begin
      function mul!(
          C::AbstractMatrix{T}, adjA::$Adj{T,<:AbstractMatrix{T}}, adjB::$Adj{T,LowerTrapezoidal{T}}) where T<:BlasFloat
        A = parent(adjA)
        B = parent(adjB)
        m, n = size(B)
        $g(view(C,:,1:n), A)
        BLAS.trmm!(
          'R', 'L', $trans, 'N', one(T), view(B.data,1:n,:), view(C,:,1:n))
        mul!(view(C,:,n+1:m), $Adj(A), $Adj(view(B.data,n+1:m,:)))
        C
      end
    end
  end

  ## UpperTrapezoidal left-multiplication

  function mul!(
      y::AbstractVector{T}, A::UpperTrapezoidal{T}, x::AbstractVector{T}) where T<:BlasFloat
    m, n = size(A)
    copyto!(y, view(x,1:m))
    lmul!(UpperTriangular(view(A.data,:,1:m)), y)
    BLAS.gemv!('N', one(T), view(A.data,:,m+1:n), view(x,m+1:n), one(T), y)
  end
  function mul!(
      C::AbstractMatrix{T}, A::UpperTrapezoidal{T}, B::AbstractMatrix{T}) where T<:BlasFloat
    m, n = size(A)
    copyto!(C, view(B,1:m,:))
    lmul!(UpperTriangular(view(A.data,:,1:m)), C)
    BLAS.gemm!(
      'N', 'N', one(T), view(A.data,:,m+1:n), view(B,m+1:n,:), one(T), C)
  end

  for (Adj, g, trans) in ((:Adjoint, :adjoint!, 'C'),
                        (:Transpose, :transpose!,  'T'))
    @eval begin
      function mul!(C::AbstractMatrix{T}, A::UpperTrapezoidal{T}, adjB::$Adj{T,<:AbstractMatrix{T}}) where T<:BlasFloat
        B = parent(adjB)
        m, n = size(A)
        $g(C, view(B,:,1:m))
        lmul!(UpperTriangular(view(A.data,:,1:m)), C)
        BLAS.gemm!(
          'N', $trans, one(T), view(A.data,:,m+1:n), view(B,:,m+1:n), one(T), C)
      end
    end
  end

  for (Adj, trans) in ((:Adjoint, 'C'), (:Transpose, 'T'))
    @eval begin
      function mul!(y::AbstractVector{T}, adjA::$Adj{T,UpperTrapezoidal{T}}, x::AbstractVector{T}) where T<:BlasFloat
        A = parent(adjA)
        m, n = size(A)
        y[1:m] = x
        BLAS.trmv!('U', $trans, 'N', view(A.data,:,1:m), view(y,1:m))
        mul!(view(y,m+1:n), $Adj(view(A.data,:,m+1:n)), x)
        y
      end
      function mul!(C::AbstractMatrix{T}, adjA::$Adj{T,UpperTrapezoidal{T}}, B::AbstractMatrix{T}) where T<:BlasFloat
        A = parent(adjA)
        m, n = size(A)
        C[1:m,:] = B
        BLAS.trmm!(
          'L', 'U', $trans, 'N', one(T), view(A.data,:,1:m), view(C,1:m,:))
        mul!(view(C,m+1:n,:), $Adj(view(A.data,:,m+1:n)), B)
        C
      end
    end
  end

  for (Adj, g, trans) in ((:Adjoint, :adjoint!, 'C'),
                        (:Transpose, :transpose!,  'T'))
    @eval begin
      function mul!(C::AbstractMatrix{T}, adjA::$Adj{T,UpperTrapezoidal{T}}, adjB::$Adj{T,AbstractMatrix{T}}) where T<:BlasFloat
        A = parent(adjA)
        B = parent(adjB)
        m, n = size(A)
        $g(view(C,1:m,:), B)
        BLAS.trmm!(
          'L', 'U', $trans, 'N', one(T), view(A.data,:,1:m), view(C,1:m,:))
        mul!(view(C,m+1:n,:), $Adj(view(A.data,:,m+1:n)), $Adj(B))
        C
      end
    end
  end

  ## UpperTrapezoidal right-multiplication

  function mul!(
      C::AbstractMatrix{T}, A::AbstractMatrix{T}, B::UpperTrapezoidal{T}) where T
    m, n = size(B)
    C[:,1:m] = A
    rmul!(view(C,:,1:m), UpperTriangular(view(B.data,:,1:m)))
    mul!(view(C,:,m+1:n), A, view(B.data,:,m+1:n))
    C
  end

  for (Adj, trans) in ((:Adjoint, 'C'), (:Transpose, 'T'))
    @eval begin
      function mul!(
          C::AbstractMatrix{T}, A::AbstractMatrix{T}, adjB::$Adj{T,UpperTrapezoidal{T}}) where T<:BlasFloat
        B = parent(adjB)
        m, n = size(B)
        copyto!(C, view(A,:,1:m))
        BLAS.trmm!('R', 'U', $trans, 'N', one(T), view(B.data,:,1:m), C)
        BLAS.gemm!(
          'N', $trans, one(T), view(A,:,m+1:n), view(B.data,:,m+1:n), one(T), C)
      end
    end
  end

  for (Adj, g) in ((:Adjoint, :adjoint!), (:Transpose, :transpose!))
    @eval begin
      function mul!(
          C::AbstractMatrix{T}, adjA::$Adj{T,<:AbstractMatrix{T}}, B::UpperTrapezoidal{T}) where T
        A = parent(adjA)
        m, n = size(B)
        $g(view(C,:,1:m), A)
        rmul!(view(C,:,1:m), UpperTriangular(view(B.data,:,1:m)))
        mul!(view(C,:,m+1:n), $Adj(A), view(B.data,:,m+1:n))
        C
      end
    end
  end

  for (Adj, g, trans) in ((:Adjoint, :adjoint!, 'C'),
                        (:Transpose, :transpose!,  'T'))
    @eval begin
      function mul!(C::AbstractMatrix{T}, adjA::$Adj{T,<:AbstractMatrix{T}}, adjB::$Adj{T,UpperTrapezoidal{T}}) where T<:BlasFloat
        A = parent(adjA)
        B = parent(adjB)
        m, n = size(B)
        mul!(C, view(A,1:m,:))
        BLAS.trmm!('R', 'U', $trans, 'N', one(T), view(B.data,:,1:m), C)
        BLAS.gemm!(
          $trans, $trans, one(T), view(A,m+1:n,:), view(B.data,:,m+1:n),
          one(T), C)
      end
    end
  end

  # standard operations

  ## left-multiplication

  function *(A::Trapezoidal{TA}, B::AbstractVector{TB}) where {TA,TB}
    T = promote_type(TA, TB)
    AT = convert(Trapezoidal{T}, A)
    BT = (T == TB ? B : convert(Array{T}, B))
    CT = Array{T}(undef, size(A,1))
    mul!(CT, AT, BT)
  end

  function *(A::Trapezoidal{TA}, B::AbstractMatrix{TB}) where {TA,TB}
    T = promote_type(TA, TB)
    AT = convert(Trapezoidal{T}, A)
    BT = (T == TB ? B : convert(Array{T}, B))
    CT = Array{T}(undef, size(A,1), size(B,2))
    mul!(CT, AT, BT)
  end
  function *(A::AbstractMatrix{TA}, B::Trapezoidal{TB}) where {TA,TB}
    T = promote_type(TA, TB)
    AT = (T == TA ? A : convert(Array{T}, A))
    BT = convert(Trapezoidal{T}, B)
    CT = Array{T}(undef, size(A,1), size(B,2))
    mul!(CT, AT, BT)
  end

  for Adj in (:Adjoint, :Transpose)
    @eval begin
      function *(adjA::$Adj{TA,<:Trapezoidal{TA}}, B::AbstractVector{TB}) where {TA,TB}
        A = parent(adjA)
        T = promote_type(TA, TB)
        AT = convert(Trapezoidal{T}, A)
        BT = (T == TB ? B : convert(Array{T}, B))
        CT = Array{T}(undef, size(A,2))
        mul!(CT, $Adj(AT), BT)
      end
      function *(adjA::$Adj{TA,<:Trapezoidal{TA}}, B::AbstractMatrix{TB}) where {TA,TB}
        A = parent(adjA)
        T = promote_type(TA, TB)
        AT = convert(Trapezoidal{T}, A)
        BT = (T == TB ? B : convert(Array{T}, B))
        CT = Array{T}(undef, size(A,2), size(B,2))
        mul!(CT, $Adj(AT), BT)
      end
      function *(A::Trapezoidal{TA}, adjB::$Adj{TA,<:AbstractMatrix{TB}}) where {TA,TB}
        B = parent(adjB)
        T = promote_type(TA, TB)
        AT = convert(Trapezoidal{T}, A)
        BT = (T == TB ? B : convert(Array{T}, B))
        CT = Array{T}(undef, size(A,1), size(B,1))
        mul!(CT, AT, $Adj(BT))
      end
      function *(adjA::$Adj{TA,<:Trapezoidal{TA}}, adjB::$Adj{TA,<:AbstractMatrix{TB}}) where {TA,TB}
        A = parent(adjA)
        B = parent(adjB)
        T = promote_type(TA, TB)
        AT = convert(Trapezoidal{T}, A)
        BT = (T == TB ? B : convert(Array{T}, B))
        CT = Array{T}(undef, size(A,1), size(B,1))
        mul!(CT, $Adj(AT), $Adj(BT))
      end

      function *(adjA::$Adj{TA,<:AbstractMatrix{TA}}, B::Trapezoidal{TB}) where {TA,TB}
        A = parent(adjA)
        T = promote_type(TA, TB)
        AT = (T == TA ? A : convert(Array{T}, A))
        BT = convert(Trapezoidal{T}, B)
        CT = Array{T}(undef, size(A,2), size(B,2))
        mul!(CT, $Adj(AT), BT)
      end
      function *(A::AbstractMatrix{TA}, adjB::$Adj{TA,<:Trapezoidal{TB}}) where {TA,TB}
        B = parent(adjB)
        T = promote_type(TA, TB)
        AT = (T == TA ? A : convert(Array{T}, A))
        BT = convert(Trapezoidal{T}, B)
        CT = Array{T}(undef, size(A,1), size(B,1))
        mul!(CT, AT, $Adj(BT))
      end
      function *(adjA::$Adj{TA,<:AbstractMatrix{TA}}, adjB::$Adj{TA,<:Trapezoidal{TB}}) where {TA,TB}
        A = parent(adjA)
        B = parent(adjB)
        T = promote_type(TA, TB)
        AT = (T == TA ? A : convert(Array{T}, A))
        BT = convert(Trapezoidal{T}, B)
        CT = Array{T}(undef, size(A,2), size(B,2))
        mul!(CT, $Adj(AT), $Adj(BT))
      end
    end
  end


text:
# load init commands
using InteractiveControlPlots, ControlSystems, GtkInteract

# create sliders etc
x = linspace(-1,1,2)
ws = logspace(-2,2,1000)
k=slider(1:0.1:3, label="K")
logT = slider(-1:0.01:1, label="T")
f1 = checkbox(false,label="Bode plot  ")
f2 = checkbox(false,label="Nyquist plot  ")
f3 = checkbox(false,label="PZ map  ")
f4 = checkbox(false,label="Step plot  ")

# create window
cg1 = cairographic(width=200,height=150)
cg2 = cairographic(width=200,height=150)
cg3 = cairographic(width=200,height=150)
cg4 = cairographic(width=200,height=150)
lT = label("T: ")

vbl(w, len) = vbox(label(w.label),hbox(hskip(len),w))
vbl(w) = vbox(label(w.label),w)

w = window(vbox(hbox(vbl(k),
    vbox(lT,logT),
    vbox(vbl(f1,25),vbl(f3,15)),
    vbox(vbl(f2,35),vbl(f4,22))),
    hbox(grow(cg1), grow(cg2)),
    hbox(grow(cg3), grow(cg4))),
    title="Example")

# loop over parameters
Reactive.preserve(map(k,logT,f1,f2,f3,f4) do k,logT,f1,f2,f3,f4
  T = 10^logT

  #Create the 2nd order system
  sys = tf(k*[1,1],[1, T, 1])

  #Create the plots
  fig1 = bode_immerse(sys, ws, 1e-2, 1e2) #Ylims = (1e-2, 1e2)
  fig2 = nyquist_immerse(sys, ws, -2, 2, -2, 2) #x and ylims
  fig3 = pzmap_immerse(sys,-4,4,-4,4) # x and ylims
  y,t,x = step(sys, 10)
  fig4 = timeplot_immerse(t, y, 0, 4) #ylims

  #Update the text
  push!(lT, "T: " * string(round(T,2)))

  #Push figures
  push!(cg1, f1 ? fig1 : nothing)
  push!(cg2, f2 ? fig2 : nothing)
  push!(cg3, f3 ? fig3 : nothing)
  push!(cg4, f4 ? fig4 : nothing)

end)

# output window for display
w


text:
using DeferredFutures
using Distributed
using Serialization
using Test

@testset "DeferredRemoteRefs" begin
    @testset "DeferredFuture Comparison" begin
        rc = RemoteChannel()
        @test DeferredFuture(rc) == DeferredFuture(rc)
        @test hash(DeferredFuture(rc)) == hash(DeferredFuture(rc))
    end

    @testset "DeferredChannel Comparison" begin
        rc = RemoteChannel()
        func = () -> RemoteChannel()
        @test DeferredChannel(rc, func) == DeferredChannel(rc, func)
        @test hash(DeferredChannel(rc, func)) == hash(DeferredChannel(rc, func))
    end

    @testset "Finalizing" begin
        df = DeferredFuture()
        finalize(df)
        @test_throws Exception isready(df)
        @test_throws Exception fetch(df)
        @test_throws Exception df[]
        @test_throws Exception put!(df, 1)
        @test_throws Exception take!(df)
        @test_throws Exception wait(df)
        finalize(df)

        dc = DeferredChannel()
        finalize(dc)
        @test_throws Exception isready(dc)
        @test_throws Exception fetch(dc)
        @test_throws Exception close(dc)
        @test_throws Exception dc[]
        @test_throws Exception put!(dc, 1)
        @test_throws Exception take!(dc)
        @test_throws Exception wait(dc)
        finalize(dc)

        dc = DeferredChannel()
        close(dc)
        @test !isready(dc)
        @test_throws Exception fetch(dc)
        @test_throws Exception dc[]
        @test_throws Exception put!(dc, 1)
        @test_throws Exception take!(dc)
        @test_throws Exception wait(dc)
        close(dc)
        finalize(dc)
        @test_throws Exception isready(dc)
        @test_throws Exception fetch(dc)
        @test_throws Exception close(dc)
        @test_throws Exception dc[]
        @test_throws Exception put!(dc, 1)
        @test_throws Exception take!(dc)
        @test_throws Exception wait(dc)

        df = DeferredFuture()
        put!(df, 1)
        @test df[] == 1
        finalize(df)
        @test_throws Exception isready(df)
        @test_throws Exception fetch(df)
        @test_throws Exception df[]
        @test_throws Exception put!(df, 1)
        @test_throws Exception take!(df)
        @test_throws Exception wait(df)

        dc = DeferredChannel()
        put!(dc, 1)
        @test dc[] == 1
        finalize(dc)
        @test_throws Exception isready(dc)
        @test_throws Exception fetch(dc)
        @test_throws Exception close(dc)
        @test_throws Exception dc[]
        @test_throws Exception put!(dc, 1)
        @test_throws Exception take!(dc)
        @test_throws Exception wait(dc)

        dc = DeferredChannel()
        put!(dc, 1)
        @test dc[] == 1
        close(dc)
        @test isready(dc)
        @test fetch(dc) == 1
        @test dc[] == 1
        @test_throws Exception put!(dc, 1)
        @test take!(dc) == 1
        @test !isready(dc)
        @test_throws Exception fetch(dc)
        @test_throws Exception dc[]
        @test_throws Exception put!(dc, 1)
        @test_throws Exception take!(dc)
        @test_throws Exception wait(dc)
        finalize(dc)
        @test_throws Exception isready(dc)
        @test_throws Exception fetch(dc)
        @test_throws Exception close(dc)
        @test_throws Exception dc[]
        @test_throws Exception put!(dc, 1)
        @test_throws Exception take!(dc)
        @test_throws Exception wait(dc)
    end

    @testset "Distributed DeferredFuture" begin
        top = myid()
        bottom = addprocs(1)[1]
        @everywhere using DeferredFutures

        try
            val = "hello"
            df = DeferredFuture(top)

            @test !isready(df)

            fut = remotecall_wait(bottom, df) do dfr
                put!(dfr, val)
            end
            @test fetch(fut) == df
            @test isready(df)
            @test fetch(df) == val
            @test wait(df) == df
            @test_throws ErrorException put!(df, val)

            @test df[] == val
            @test df[5] == 'o'

            @test df.outer.where == top
            @test fetch(df.outer).where == bottom

            reset!(df)
            @test !isready(df)
            put!(df, "world")
            @test fetch(df) == "world"
            finalize(df)
            @test_throws Exception isready(df)
            @test_throws Exception fetch(df)
            @test_throws Exception df[]
            @test_throws Exception put!(df, 1)
            @test_throws Exception take!(df)
            @test_throws Exception wait(df)
            finalize(df)
        finally
            rmprocs(bottom)
        end
    end

    @testset "Distributed DeferredChannel" begin
        top = myid()
        bottom = addprocs(1)[1]
        @everywhere using DeferredFutures

        try
            val = "hello"
            channel = DeferredChannel(top, 32)

            @test !isready(channel)

            fut = remotecall_wait(bottom, channel) do dfr
                put!(dfr, val)
            end
            @test fetch(fut) == channel
            @test isready(channel)
            @test fetch(channel) == val
            @test wait(channel) == channel

            @test channel[] == val
            @test channel[5] == 'o'

            put!(channel, "world")
            @test take!(channel) == val
            @test fetch(channel) == "world"

            @test channel.outer.where == top
            @test fetch(channel.outer).where == bottom

            reset!(channel)
            @test !isready(channel)
            put!(channel, "world")
            @test fetch(channel) == "world"
            finalize(channel)
            @test_throws Exception isready(channel)
            @test_throws Exception fetch(channel)
            @test_throws Exception close(channel)
            @test_throws Exception channel[]
            @test_throws Exception put!(channel, 1)
            @test_throws Exception take!(channel)
            @test_throws Exception wait(channel)
            finalize(channel)
        finally
            rmprocs(bottom)
        end
    end

    @testset "Allocation" begin
        rand_size = 800000000  # sizeof(rand(10000, 10000))
        GC.gc()
        main_size = Base.summarysize(Distributed)

        top = myid()
        bottom = addprocs(1)[1]
        @everywhere using DeferredFutures

        try
            df = DeferredFuture(top)

            remote_size = remotecall_fetch(bottom, df) do dfr
                GC.gc()
                main_size = Base.summarysize(Distributed)

                # the DeferredFuture is initialized and the data is stored on bottom
                put!(dfr, rand(10000, 10000))
                main_size
            end

            GC.gc()
            # tests that the data has not been transfered to top
            @test Base.summarysize(Distributed) < main_size + rand_size

            remote_size_new = remotecall_fetch(bottom) do
                GC.gc()
                Base.summarysize(Distributed)
            end

            # tests that the data still exists on bottom
            @test remote_size_new >= remote_size + rand_size
        finally
            rmprocs(bottom)
        end
    end

    @testset "Transfer" begin
        rand_size = 800000000  # sizeof(rand(10000, 10000))
        GC.gc()
        main_size = Base.summarysize(Main)

        top = myid()
        left, right = addprocs(2)
        @everywhere using DeferredFutures

        try
            df = DeferredFuture(top)

            left_remote_size = remotecall_fetch(left, df) do dfr
                GC.gc()
                main_size = Base.summarysize(Main)
                put!(dfr, rand(10000, 10000))
                main_size
            end

            right_remote_size = remotecall_fetch(right, df) do dfr
                GC.gc()
                main_size = Base.summarysize(Main)
                global data = fetch(dfr)
                main_size
            end

            GC.gc()
            @test Base.summarysize(Main) < main_size + rand_size

            right_remote_size_new = remotecall_fetch(right) do
                GC.gc()
                Base.summarysize(Main)
            end

            @test right_remote_size_new >= right_remote_size + rand_size
        finally
            rmprocs([left, right])
        end
    end

    @testset "@defer" begin
        ex = macroexpand(@__MODULE__, :(@defer RemoteChannel(()->Channel(5))))
        ex = macroexpand(@__MODULE__, :(@defer RemoteChannel()))

        channel = @defer RemoteChannel(()->Channel(32))

        put!(channel, 1)
        put!(channel, 2)

        @test fetch(channel) == 1
        @test take!(channel) == 1
        @test fetch(channel) == 2

        fut = macroexpand(@__MODULE__, :(@defer Future()))
        other_future = macroexpand(@__MODULE__, :(@defer Future()))

        @test_throws LoadError macroexpand(@__MODULE__, :(@defer mutable struct Foo end))
        try
            macroexpand(@__MODULE__, :(@defer mutable struct Foo end))
            @test false
        catch e
            @test e.error isa AssertionError
        end

        @test_throws LoadError macroexpand(@__MODULE__, :(@defer Channel()))
        try
            macroexpand(@__MODULE__, :(@defer Channel()))
            @test false
        catch e
            @test e.error isa AssertionError
        end
        close(channel)
    end

    @testset "Show" begin
        rc = RemoteChannel()
        rc_params = "($(rc.where),$(rc.whence),$(rc.id))"

        @test sprint(show, DeferredFuture(rc)) == "DeferredFuture at $rc_params"

        dc = DeferredChannel(rc, print)
        @test sprint(show, dc) == "DeferredChannel(print) at $rc_params"
    end

    @testset "Serialization" begin
        @testset "DeferredFuture serialization on same process" begin
            df = DeferredFuture(myid())

            io = IOBuffer()
            serialize(io, df)
            seekstart(io)
            deserialized_df = deserialize(io)
            close(io)

            @test deserialized_df == df
        end

        @testset "DeferredFuture serialization on a cluster" begin
            df1 = DeferredFuture(myid())
            df2 = DeferredFuture(myid())

            io = IOBuffer()
            serialize(io, df1)
            df1_string = take!(io)
            close(io)

            put!(df2, 28)

            io = IOBuffer()
            serialize(io, df2)
            df2_string = take!(io)
            close(io)

            bottom = addprocs(1)[1]
            @everywhere using DeferredFutures
            @everywhere using Serialization

            df3_string = ""
            try
                df3_string = @fetchfrom bottom begin
                    io = IOBuffer()
                    write(io, df1_string)
                    seekstart(io)
                    bottom_df1 = deserialize(io)
                    close(io)

                    put!(bottom_df1, 37)

                    io = IOBuffer()
                    write(io, df2_string)
                    seekstart(io)
                    bottom_df2 = deserialize(io)
                    close(io)

                    @test isready(bottom_df2) == true
                    @test fetch(bottom_df2) == 28
                    reset!(bottom_df2)

                    df3 = DeferredFuture(myid())
                    put!(df3, 14)

                    io = IOBuffer()
                    serialize(io, df3)
                    df3_string = take!(io)
                    close(io)

                    return df3_string
                end

                @test isready(df1) == true
                @test fetch(df1) == 37

                @test isready(df2) == false

                @test df3_string != ""

            finally
                rmprocs(bottom)
            end

            io = IOBuffer()
            write(io, df3_string)
            seekstart(io)
            bottom_df3 = deserialize(io)
            close(io)

            @test_broken isready(bottom_df3) == true
            @test_broken fetch(bottom_df3) == 14
        end

        @testset "DeferredChannel serialization on same process" begin
            dc = DeferredChannel()

            io = IOBuffer()
            serialize(io, dc)
            seekstart(io)
            deserialized_dc = deserialize(io)
            close(io)

            @test deserialized_dc == dc
        end

        @testset "DeferredChannel serialization on a cluster" begin
            dc1 = DeferredChannel()
            dc2 = DeferredChannel()

            io = IOBuffer()
            serialize(io, dc1)
            dc1_string = take!(io)
            close(io)

            put!(dc2, 28)

            io = IOBuffer()
            serialize(io, dc2)
            dc2_string = take!(io)
            close(io)

            bottom = addprocs(1)[1]
            @everywhere using DeferredFutures
            @everywhere using Serialization

            dc3_string = ""
            try
                dc3_string = @fetchfrom bottom begin
                    io = IOBuffer()
                    write(io, dc1_string)
                    seekstart(io)
                    bottom_dc1 = deserialize(io)
                    close(io)

                    put!(bottom_dc1, 37)

                    io = IOBuffer()
                    write(io, dc2_string)
                    seekstart(io)
                    bottom_dc2 = deserialize(io)
                    close(io)

                    @test isready(bottom_dc2) == true
                    @test fetch(bottom_dc2) == 28
                    reset!(bottom_dc2)

                    dc3 = DeferredChannel()
                    put!(dc3, 14)

                    io = IOBuffer()
                    serialize(io, dc3)
                    dc3_string = take!(io)
                    close(io)

                    return dc3_string
                end

                @test isready(dc1) == true
                @test fetch(dc1) == 37

                @test isready(dc2) == false

                @test dc3_string != ""

            finally
                rmprocs(bottom)
            end

            io = IOBuffer()
            write(io, dc3_string)
            seekstart(io)
            bottom_dc3 = deserialize(io)
            close(io)

            @test_broken isready(bottom_dc3) == true
            @test_broken fetch(bottom_dc3) == 14
        end

        @testset "DeferredFuture serialization as part of another object" begin
            pnums = addprocs(1)
            @everywhere using DeferredFutures

            try
                x = ()->3
                df = DeferredFuture()
                result = @fetch (df, x)

                @test result[1] == df
                @test result[2]() == 3
            finally
                rmprocs(pnums)
            end
        end

        @testset "DeferredChannel serialization as part of another object" begin
            pnums = addprocs(1)
            @everywhere using DeferredFutures

            try
                x = ()->3
                dc = DeferredChannel()
                result = @fetch (dc, x)

                @test result[1] == dc
                @test result[2]() == 3
            finally
                rmprocs(pnums)
            end
        end
    end
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 IoK8sApiAuditregistrationV1alpha1Webhook <: SwaggerModel
    clientConfig::Any # spec type: Union{ Nothing, IoK8sApiAuditregistrationV1alpha1WebhookClientConfig } # spec name: clientConfig
    throttle::Any # spec type: Union{ Nothing, IoK8sApiAuditregistrationV1alpha1WebhookThrottleConfig } # spec name: throttle

    function IoK8sApiAuditregistrationV1alpha1Webhook(;clientConfig=nothing, throttle=nothing)
        o = new()
        validate_property(IoK8sApiAuditregistrationV1alpha1Webhook, Symbol("clientConfig"), clientConfig)
        setfield!(o, Symbol("clientConfig"), clientConfig)
        validate_property(IoK8sApiAuditregistrationV1alpha1Webhook, Symbol("throttle"), throttle)
        setfield!(o, Symbol("throttle"), throttle)
        o
    end
end # type IoK8sApiAuditregistrationV1alpha1Webhook

const _property_map_IoK8sApiAuditregistrationV1alpha1Webhook = Dict{Symbol,Symbol}(Symbol("clientConfig")=>Symbol("clientConfig"), Symbol("throttle")=>Symbol("throttle"))
const _property_types_IoK8sApiAuditregistrationV1alpha1Webhook = Dict{Symbol,String}(Symbol("clientConfig")=>"IoK8sApiAuditregistrationV1alpha1WebhookClientConfig", Symbol("throttle")=>"IoK8sApiAuditregistrationV1alpha1WebhookThrottleConfig")
Base.propertynames(::Type{ IoK8sApiAuditregistrationV1alpha1Webhook }) = collect(keys(_property_map_IoK8sApiAuditregistrationV1alpha1Webhook))
Swagger.property_type(::Type{ IoK8sApiAuditregistrationV1alpha1Webhook }, name::Symbol) = Union{Nothing,eval(Meta.parse(_property_types_IoK8sApiAuditregistrationV1alpha1Webhook[name]))}
Swagger.field_name(::Type{ IoK8sApiAuditregistrationV1alpha1Webhook }, property_name::Symbol) =  _property_map_IoK8sApiAuditregistrationV1alpha1Webhook[property_name]

function check_required(o::IoK8sApiAuditregistrationV1alpha1Webhook)
    (getproperty(o, Symbol("clientConfig")) === nothing) && (return false)
    true
end

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


text:
_axisdim(a, ax) = axisdim(a, ax)
_axisdim(a, ax::Integer) = ax

function _makewin(fid, axis, windowfunc)
    dim = _axisdim(fid, axis)
    dimlen = size(fid, dim)
    w = windowfunc.((0:dimlen-1)/dimlen)
    # Reshape to broadcast `w` only along time dimension
    windowshape = ntuple(i->i==dim ? dimlen : 1, ndims(fid))
    reshape(w, windowshape...)
end

"""
Apply `windowfunc` over the dimensionless time range `(0:tlen-1)/tlen)` of
`fid` along dimension `axis` which can be an `Axis` or `Integer`.

See also `apply_window`.
"""
function apply_window!(fid::AxisArray, axis, windowfunc)
    w = _makewin(fid, axis, windowfunc)
    fid .*= w
    fid
end

"""
    apply_window(signal, axis1=>win1, axis2=>win2, ...)

Apply window functions `win1` along `axis1`, `win2` along `axis2`, etc.  The
windows should be functions over the dimensionless time domain `0..1`; the axes
should be of type `AxisArrays.Axis`.

# Background

For the purposes of spectroscopy, the MR signal is generally an oscillating
decaying signal in the time domain, with different metabolites having different
decay constants. When computing a spectrum with a simple Fourier transform, one
may window the time domain data to emphasize one metabolite or another in the
resulting spectrum. (Spatially resolved spectroscopy involves some kind of MR
echo so the initial part of the signal *may* be increasing, but the this doesn't
affect the irreversible decay mechanisms.)

For *qualitative* exploration of the spectra, one may therefore want to choose
out of several windowing functions with convenient parametric forms. Note that
for *quantitative* measurement of metabolite concentration, windowing distorts
the relative peak volumes and it's probably better to choose a time domain
fitting method.

Note that these MR-specific reasons are subtly different - and apply in addition
to - the generic signal processing reasons for time domain windowing.
"""
function apply_window(fid::AxisArray, windows::Pair...)
    fid = copy(fid)
    axis,winfunc = windows[1]
    for (axis,winfunc) in windows
        apply_window!(fid, axis, winfunc)
    end
    fid
end

"""
    sinebell(t; skew=1, pow=1)

The sinebell window on the dimensionless time range `t ∈ 0..1`.

Use `pow=2` for sinebell squared. Setting the `skew` parameter to something
other than the default of `1.0` gives a skewed sinebell window, as in Felix NMR.
"""
function sinebell(t::Number; skew=1.0, pow=1.0)
    # Guess at the Felix parameterization via inspection of the window in the
    # UI; seen via "Open and Process"->"Window"->"Real time"
    sin(pi*t^skew)^pow
end

max_stars_repo_name: TRIImaging/MagneticResonanceSignals

max_issues_repo_name: TRIImaging/MagneticResonanceSignals

max_forks_repo_name: TRIImaging/MagneticResonanceSignals


text:
module effective_interaction_simple

using ..system, ..interaction

const origin = zeros(Float64, 3)
const e_z = Float64[0.,0.,1.]
const gamma = 1.

function effective_interactions(S::SpinCollection)
    omega_eff::Float64 = 0.
    gamma_eff::Float64 = 0.
    for spin = S.spins
        omega_eff += interaction.Omega(origin, spin.position, S.polarizations[1], S.polarizations[1], gamma, gamma)
        gamma_eff += interaction.Gamma(origin, spin.position, S.polarizations[1], S.polarizations[1], gamma, gamma)
    end
    return omega_eff, gamma_eff
end

# Finite symmetric systems

function triangle_orthogonal(a::Float64)
    positions = Vector{Float64}[[a,0,0], [a/2, sqrt(3. /4)*a,0]]
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

function square_orthogonal(a::Float64)
    positions = Vector{Float64}[[a,0,0], [0,a,0], [a,a,0]]
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

function polygon_orthogonal(N::Int, a::Float64)
    @assert N>2
    dα = 2*pi/N
    R = a/(2*sin(dα/2))
    positions = Vector{Float64}[]
    for i=1:(N-1)
        x = R*cos(i*dα)
        y = R*sin(i*dα)
        push!(positions, Float64[x-R,y,0])
    end
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

function rectangle_orthogonal(a::Float64, b::Float64)
    positions = Vector{Float64}[[a,0,0], [0,b,0], [a,b,0]]
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

function cube_orthogonal(a::Float64)
    positions = Vector{Float64}[]
    for ix=0:1, iy=0:1, iz=0:1
        if ix==0 && iy==0 && iz==0
            continue
        end
        push!(positions, [ix*a, iy*a, iz*a])
    end
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

function box_orthogonal(a::Float64, b::Float64, c::Float64)
   positions = Vector{Float64}[]
    for ix=0:1, iy=0:1, iz=0:1
        if ix==0 && iy==0 && iz==0
            continue
        end
        push!(positions, [ix*a, iy*b, iz*c])
    end
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

# Infinite 1D symmetric systems

function chain(a::Float64, Θ, N::Int)
    positions = Vector{Float64}[]
    for ix=-N:N
        if ix==0
            continue
        end
        push!(positions, [ix*a, 0., 0.])
    end
    S = SpinCollection(positions, Float64[cos(Θ), 0., sin(Θ)]; gammas=gamma)
    return effective_interactions(S)
end

function chain_orthogonal(a::Float64, N::Int)
    positions = Vector{Float64}[]
    for ix=-N:N
        if ix==0
            continue
        end
        push!(positions, [ix*a, 0., 0.])
    end
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

# Infinite 2D symmetric systems

function squarelattice_orthogonal(a::Float64, N::Int)
    positions = Vector{Float64}[]
    for ix=-N:N, iy=-N:N
        if ix==0 && iy==0
            continue
        end
        push!(positions, [ix*a, iy*a, 0.])
    end
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

function hexagonallattice_orthogonal(a::Float64, N::Int)
    positions = Vector{Float64}[]
    ax = sqrt(3.0/4)*a
    for iy=1:N
        push!(positions, [0, iy*a, 0])
        push!(positions, [0, -iy*a, 0])
    end
    for ix=[-1:-2:-N; 1:2:N]
        Ny = div(2*N+1-abs(ix),2)
        for iy=0:Ny-1
            push!(positions, [ax*ix, (0.5+iy)*a, 0])
            push!(positions, [ax*ix, -(0.5+iy)*a, 0])
        end
    end
    for ix=[-2:-2:-N; 2:2:N]
        Ny = div(2*N-abs(ix),2)
        push!(positions, [ax*ix, 0, 0])
        for iy=1:Ny
            push!(positions, [ax*ix, iy*a, 0])
            push!(positions, [ax*ix, -iy*a, 0])
        end
    end
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

# Infinite 3D symmetric systems

function cubiclattice_orthogonal(a::Float64, N::Int)
    positions = Vector{Float64}[]
    for ix=-N:N, iy=-N:N, iz=-N:N
        if ix==0 && iy==0 && iz==0
            continue
        end
        push!(positions, [ix*a, iy*a, iz*a])
    end
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

function tetragonallattice_orthogonal(a::Float64, b::Float64, N::Int)
    positions = Vector{Float64}[]
    for ix=-N:N, iy=-N:N, iz=-N:N
        if ix==0 && iy==0 && iz==0
            continue
        end
        push!(positions, [ix*a, iy*a, iz*b])
    end
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

function hexagonallattice3d_orthogonal(a::Float64, b::Float64, N::Int)
    positions = Vector{Float64}[]
    ax = sqrt(3.0/4)*a
    for iz=-N:N
        for iy=1:N
            push!(positions, [0, iy*a, iz*b])
            push!(positions, [0, -iy*a, iz*b])
        end
        for ix=[-1:-2:-N; 1:2:N]
            Ny = div(2*N+1-abs(ix),2)
            for iy=0:Ny-1
                push!(positions, [ax*ix, (0.5+iy)*a, iz*b])
                push!(positions, [ax*ix, -(0.5+iy)*a, iz*b])
            end
        end
        for ix=[-2:-2:-N; 2:2:N]
            Ny = div(2*N-abs(ix),2)
            push!(positions, [ax*ix, 0, iz*b])
            for iy=1:Ny
                push!(positions, [ax*ix, iy*a, iz*b])
                push!(positions, [ax*ix, -iy*a, iz*b])
            end
        end
    end
    S = SpinCollection(positions, e_z; gammas=gamma)
    return effective_interactions(S)
end

end # module


text:
cd(@__DIR__)
using Plots, Flux, LinearAlgebra, Parameters, Statistics, Random, Printf, OrdinaryDiffEq, DSP, SparseArrays, Dates
using Flux: params, jacobian, train!, data
using Serialization: serialize, deserialize
default(lab="", grid=false)
# include("SkipRNN.jl")

@userplot Eigvalplot
@recipe function eigvalplot(A::Eigvalplot)
    e = data(A.args[1]) |> eigvals# .|> log
    title --> "Eigenvalues"
    ar = range(0, stop=2pi, length=40)
    @series begin
        linestyle := :dash
        color := :black
        cos.(ar), sin.(ar)
    end
    @series begin
        seriestype := :scatter
        real.(e), imag.(e)
    end
end

const h = 0.1

function pendcart(xd,x,p,t)
    g = 9.82; l = 1.0; d = 0.5
    u = 0
    xd[1] = x[2]
    xd[2] = -g/l * sin(x[1]) + u/l * cos(x[1]) - d*x[2]
    xd
end

function generate_data(T=20)
    u0 = @. [pi, 6] * 2 *(rand()-0.5)
    tspan = (0.0,T)
    prob = ODEProblem(pendcart,u0,tspan)
    sol = solve(prob,Tsit5())
    z = reduce(hcat, sol(0:h:T).u)
    y = vcat(sin.(z[1:1,:]), cos.(z[1:1,:])) .+ 0.05 .* randn.()
    z,y
end

##
# function dostuff()
trajs = [generate_data(10) for i = 1:50]
const dataset = [(t[2], t[1]) for t in trajs]

const nz = 3
const nx = size(trajs[1][2], 1)
const nu = 0
const nα = 5

const pxz = Chain(Dense(nz,30,relu), Dense(30,nx))
const f = Chain(Dense(nx+nz+nu,30,relu), Dense(30,2*nz))
const αnet = Chain(Dense(nz,nα), softmax)
# q = Chain(Dense())
const initialw = Chain(Dense(3nx,30,relu), Dense(30,2*nz))
const initialz = Chain(Dense(nz,30,relu), Dense(30,nz))

const A = [param(0.001randn(nz,nz)) for _ = 1:nα]
const C = [param(0.001randn(nz,nz)) for _ = 1:nα]
const lσp = zeros(nz) # log space
const μp  = zeros(nz)
const lσx = param(zeros(nx)) # log space

pars = params((
pxz,
f,
αnet,
initialw,
initialz,
(A...),
(C...),
lσx,
))

##
const c = [0.2]

function finitialw(x)
    μσw = initialw(x[:,1:3][:])
    μw  = μσw[1:end÷2]
    σw = μσw[(end÷2+1):end]
    varw = σw .^ 2 .+ 1e-5
    w   = μw .+ σw .* randn(nz)

    return w, μw, varw
end

function lossindividual(x,y)
    w, μw, σw² = finitialw(x)
    z  = Float64.(initialz(w))
    l0 = kl_normalgauss(μw, σw², c[])
    σx  = exp.(lσx) .+ 1e-5
    sum(1:size(x,2)-1) do i
        xi  = x[:,i]
        # xi1 = x[:,i+1]
        μx  = pxz(z)
        lrec  = sum(abs2(xi[i] - μx[i])/(2.0 * σx[i]^2) + 0.5 * log(2π) + lσx[i]  for i in eachindex(xi))
        μσq = f([z; TrackedArray(xi)]) # TODO: I change from x[t+1] to x[t]
        μq  = μσq[1:end÷2]
        σq  = μσq[(end÷2+1):end]
        w   = μq .+ σq .* randn(nz)

        α  = αnet(z)
        Ai = sum(α[i]*A[i] for i = 1:nα)
        Ci = sum(α[i]*C[i] for i = 1:nα)
        z  = Ai*z + Ci*w
        lkl = kl_normalgauss(μq, σq .+ 1e-5, c[])
         [c[] * lrec + lkl, lrec, kl_normalgauss(μq, σq .+ 1e-5)]
    end + [l0, 0.0, l0]
end

function loss(x,y)
    lossindividual(x,y)[1]
end

function kl_normalgauss(μ1, σ1, c = 1)
    l2π = log(2π)
    lσ1 = log.(σ1)
    0.5 * sum(c*(l2π) - (l2π + 2lσ1[i]) +
    c*(exp(2lσ1[i]) + abs2(μ1[i])) - 1.0 for i in eachindex(μ1))
end

function kl(μ1, lσ1, μ2, lσ2, c = 1)
    l2π = log(2π)
    0.5 * sum(c*(l2π + 2lσ2[i]) - (l2π + 2lσ1[i]) +
    c*(exp(2lσ1[i]) + abs2(μ1[i] - μ2[i]))/(exp(2lσ2[i]) + 1e-5) - 1.0 for i in eachindex(μ1))
end

simulate(T,x::AbstractMatrix) = simulate(T, finitialw(x)[1])

function simulate(T, w = exp.(lσp) .* randn(nz))
    z,x = [],[]
    z = [data(initialz(w))]
    for i = 1:T
        μx = pxz(z[end])
        push!(x, data(μx))# .+ exp.(σx) .* randn(nx)))

        μσq = f([z[end]; x[end]]) # Wrong time index of xsample?
        μq = μσq[1:end÷2]
        σq = μσq[(end÷2+1):end]
        wsample = μq# .+ σq .* randn(nz)

        α = αnet(z[end])
        Ai = sum(α[i]*A[i] for i = 1:nα)
        Ci = sum(α[i]*C[i] for i = 1:nα)
        push!(z, data(Ai*z[end] + Ci*wsample))
    end
    copy(hcat(z...)'),copy(hcat(x...)')
end

function dvbf(x)
    w = finitialw(x)[1]
    z,xf = [],[]
    z = [data(initialz(w))]
    for i = 1:size(x,2)
        μx = pxz(z[end])
        push!(xf, data(μx))# .+ exp.(σx) .* randn(nx)))
        μσq = f([z[end]; x[:,i]]) # Wrong time index of xsample?
        μq = μσq[1:end÷2]
        σq = μσq[(end÷2+1):end]
        wsample = μq# .+ σq .* randn(nz)

        α = αnet(z[end])
        Ai = sum(α[i]*A[i] for i = 1:nα)
        Ci = sum(α[i]*C[i] for i = 1:nα)
        push!(z, data(Ai*z[end] + Ci*wsample))
    end
    copy(hcat(z[1:end-1]...)'),copy(hcat(xf...)')

end

losses = [lossindividual(dataset[1]...)]
opt = ADADelta(0.1)
# Flux.Optimise.updaterule(opt, pars)
@progress for i = 1:2000
    global c
    # initializer()
    train!(loss, pars, dataset, opt)
    Ta  = 100
    c[] = min(1, c[] + 1/Ta)

    if (i < 50 && i % 5 == 0) || i % 25 == 0
        l = mean(1:50) do i
            data(lossindividual(dataset[i]...))
        end
        push!(losses, l)
        @printf("Iter : %d  Loss: %8.4g Rec: %8.4g KL: %8.4g\n", i, data(l[1]), data(l[2]), data(l[3]));
        #plot(losses, layout=1, subplot=1, xscale=:log10, size=(600,400)) |> display
        #     # plot!(dataset[1][2]', subplot=2, c=:blue)
        #     # plot!(sim(dataset[1][1])', subplot=2, c=:orange)
        #     # plot!(data(model(dataset[1][1])'), subplot=2, c=:green, ylims=extrema(dataset[1][2]))
        #     #
        #     # # plot!(data(model(dataset[1][1])'), subplot=1, lab="yhat")
        #     # e = data(rnn.cell.Wh) |> eigvals# .|> log
        #     # scatter!(real.(eA), imag.(eA), lab="Eigenvalues A", subplot=3)
        #     # scatter!(real.(e), imag.(e), lab="Eigenvalues Wh", subplot=3)
        #     # alpha = range(0, stop=2pi, length=30)
        #     # plot!(cos.(alpha), sin.(alpha), l=:dash, subplot=3)#|> display
        #     # gui()
        fo = open("/local/work/fredrikb/res_$(splitdir(@__FILE__)[2])_$(now())", "w");
        serialize(fo, (trajs, dataset, pxz, f, αnet, initialw, initialz, A, C, lσx, opt));
        close(fo);
    end
end
println("Done")

fo = open("/local/work/fredrikb/res_DVBF.jl_2018-10-21T14:21:40.489")
( trajs, dataset, pxz, f, αnet, initialw, initialz, A, C, lσx, opt) = deserialize(fo)
##
z,x = simulate(300)
plot3d(mapslices(x->[x],z,dims=1)..., c=:orange, grid=true)
# plot(losses)
##
i = rand(1:length(dataset))
z,x = simulate(100, dataset[i][1])
plot(z, l=(:orange, ), lab="Filter state")
plot!(x, c=:red, lab="Filter meas")
plot!(dataset[i][2]', c=:blue, lab="State")
plot!(dataset[i][1]', c=:green, lab="Meas")

##

i = rand(1:length(dataset))
θ = dataset[i][2][1,:]
θv = dataset[i][2][2,:]
zf,xf = dvbf(dataset[i][1])

plot(trajs[i][1]', l=(:blue,), layout=2, subplot=1)
plot!(zf, l=(:dash,:orange), subplot=1)
plot!(xf, l=(:red, :dash), subplot=1)
scatter3d!(zf[:,1],zf[:,2],zf[:,3], subplot=2, zcolor=θ)

display(current())

##
plotangle = plot()
plotvelocity = plot()
for i in eachindex(dataset)
    θ = dataset[i][2][1,2:end]
    θv = dataset[i][2][2,2:end]
    zf,xf = dvbf(dataset[i][1])
    scatter3d!(plotangle, zf[2:end,1],zf[2:end,2],zf[2:end,3], zcolor=θ)
    scatter3d!(plotvelocity, zf[2:end,1],zf[2:end,2],zf[2:end,3], zcolor=θv)
end
display(current())
savefig(plotangle, "angle.pdf")
savefig(plotvelocity, "velocity.pdf")

##

θa = Float64[]; θva = Float64[]; zfa = []
for i in eachindex(dataset)
    global θa, θva, zfa
    θ = dataset[i][2][1,:]
    θv = dataset[i][2][2,:]
    zf,xf = dvbf(dataset[i][1])
    θa = vcat(θa, θ)
    θva = vcat(θva, θv)
    push!(zfa, zf)
end
zfa = reduce(vcat, zfa)

k = zfa\θa
kv = zfa\θva
scatter(θa, zfa)
scatter(zfa*k, θa, layout=2, xlabel="\$\\theta\$")
scatter!(zfa*kv, θva, subplot=2, xlabel="\$\\dot\\theta\$")

##
plot()
foreach(eachindex(A)) do i
    eigvalplot!(A[i])
end
display(current())

# avec = range(0, stop=2pi, length=50)
# vvec = range(-3, stop=3, length=50)
# Avec = repeat(avec, 1, 50)
# Vvec = repeat(vvec', 50, 1)


text:
@testset "errors" begin
    # Invalid usage of the macro
    try
        @aml function(x)
        end
    catch
    end

    # Already defined type field name
    try
        @aml mutable struct yy "~"
            n, "~"
        end

        @aml mutable struct tt "~"
            yy, "~"
        end
    catch
    end

end


text:
using BinDeps

if is_windows()
    error("The ZStd package does not currently support Windows.")
end

BinDeps.@setup

vers = v"1.1.0"

zstd = library_dependency("zstd", aliases=["libzstd", "libzstd.$(vers.major)", "libzstd.$vers"])

if is_apple()
    using Homebrew
    provides(Homebrew.HB, "zstd", zstd, os=:Darwin)
else
    provides(Sources, URI("https://github.com/facebook/zstd/archive/v$vers.tar.gz"), zstd,
             unpacked_dir="zstd-$vers")

    provides(BuildProcess, (@build_steps begin
        GetSources(zstd)
        @build_steps begin
            ChangeDirectory(joinpath(srcdir(zstd), "zstd-$vers"))
            FileRule(joinpath(libdir(zstd), "libzstd." * BinDeps.shlib_ext), @build_steps begin
                CreateDirectory(libdir(zstd))
                `make install PREFIX=$(BinDeps.usrdir(zstd))`
            end)
        end
    end), zstd)
end

BinDeps.@install Dict(:zstd => :libzstd)


text:
function get_VAR_lag_length(D::Array, pbar::Integer, ic::ASCIIString, i::Bool=true)
        IC   = zeros(pbar,1)
        for p = 1:pbar
            i==true ? V = VAR(D,p,true) : V = VAR(D,p,false)
            n,m = size(V.X)
            t = n-p
            sig = V.Σ/t
        if ic == "aic"
            IC[p] = log(det(sig))+2*p*m^2/t
        elseif ic == "bic"
            IC[p] = log(det(sig))+(m^2*p)*log(t)/t
        elseif ic == "hqc"
            IC[p] = log(det(sig))+2*log(log(t))*m^2*p/t
        elseif ic == "aicc"
            b = t/(t-(p*m+m+1))
            IC[p] = t*(log(det(sig))+m)+2*b*(m^2*p+m*(m+1)/2)
        elseif error("'ic' must be aic, bic, hqc or aicc")
        end
     end
     length_ic = indmin(IC)
     println("Using $ic the best lag-length is $length_ic")
     return length_ic
end