text:
function thinSt(X::Array{T}, layers::stParallel; nonlinear::Function=abs, subsam::Function=bspline, stType::String="full", outputSubsample::Tuple{Int,Int}=(-1,1), totalScales = [-1 for i=1:layers.m+1]) where {T<:Real}
  @assert length(totalScales)==layers.m+1
  # Insist that X has to have at least one extra meta-dimension, even if it is 1
  if length(size(X)) == length(size(layers.subsampling))
    X = reshape(X, (1,size(X)...));
  end
  numChildλ = 0
  #println("initial")
  k, n, q, dataSizes, outputSizes, resultingSize = calculateThinStSizes(layers, outputSubsample, size(X), totalScales = totalScales)
  #println("the required sizes are dataSizes = $(dataSizes)")
  #println("outputSizes = $(outputSizes)")
  nextData = [reshape(X, (size(X)..., 1)),]
  # create list of size references if we're subsampling the output an extra amount

  if outputSubsample[1] > 1 || outputSubsample[2] > 1
    concatOutput = SharedArray(zeros(Float64, outputSizes[1][1:end-k-1]..., sum(resultingSize.*q)))
  else
    concatOutput = SharedArray(zeros(Float64, outputSizes[1][1:end-k-1]..., sum([prod(x[end-k:end]) for x in outputSizes])...))
  end
  # keep track of where we are in each meta dimension
  outPos = ones(Int64, outputSizes[1][1:end-k-1]...)
  nScalesLayers = [numScales(layers.shears[i], max(dataSizes[i][1]-1,1)) for i=1:length(layers.shears)]
  for (i,layer) in enumerate(layers.shears[1:layers.m])
    #println("starting the actual transformation in layer $(i)")
    cur = nextData #data from the previous layer
    #println("size(cur) = $(size(cur[1]))")
    # only store the intermediate results in intermediate layers
    if i < layers.m
      # store the data in a channel, along with the number of it's parent
      dataChannel = Channel{Tuple{Array{length(dataSizes[i+1])},Int64}}(size(cur)[end])
      nextData = Array{SharedArray{Float64, length(dataSizes[i+1])},1}(undef, size(cur[1])[end]*length(cur))
    end
    # if we're decreasing, we need to keep track of what path the current index corresponds to
    if stType=="decreasing" && i>=2
      keeper = [1 for k=1:i-1]
    elseif stType=="decreasing" && i==1
      keeper = [1 for k=1:i]
    end

    ######################################################################################################################################################################################################################################################################################
    innerAxes = axes(cur[1])[end-k:end-1] # effectively a set of colons of length k, to be used for input
    innerSub = (Base.OneTo(x) for x in outputSizes[i][end-k:end-1]) # effectively a set of colons of length k, to be used for the subsampled output

    outerAxes = axes(cur[1])[1:(end-k-1)] # the same idea as innerAxes, but for the example indices
    # precompute the wavelets
    if stType=="decreasing" && i>=2
      numChildλ = numChildren(keeper, layers, nScalesLayers)
      daughters = computeWavelets(size(cur[1])[end-1], layers.shears[i], nScales=numChildλ)
    else
      daughters = computeWavelets(size(cur[1])[end-1], layers.shears[i],nScales=totalScales[i])
    end
    # TODO: compute a plan_fft that is efficient
    # compute the final mother if we're on the last layer
    if i==layers.m
      #println("On i=$(i), with nScales=1")
      mother = computeWavelets(size(cur[1])[end-1], layers.shears[i+1], nScales=1)
      concatStartLast = sum(q[1:end-1].*resultingSize[1:end-1])+1
      #println("last layer starting at $(concatStartLast)")
      nScalesLastLayer = getQ(layers, n, totalScales; product=false)[end-1]
      #println("nScalesLastLayer = $(nScalesLastLayer)")
    end
    concatStart = reduce(+, q[1:i-1].*resultingSize[1:i-1],init=0)+1
    #println("starting at $(concatStart)")
    listOfProcessResults = Array{Future, 1}(undef, size(cur[1])[end]*length(cur)) # parallel processing tools; the results are spread over the previous 2 layers paths
    # iterate over the scales two layers back
    for (j,x) in enumerate(cur)
      # iterate over all paths in layer i
      for λ = 1:size(x)[end]
        # do the actual transformation; the last layer requires far less saved than the mid-layers
        if i<layers.m
          if stType == "decreasing"
             futureThing = @spawn transformMidLayer!(layers, concatOutput, x, nextData, outerAxes, innerAxes, innerSub, k, i, nScalesLayers, daughters, totalScales[i], nonlinear, subsam, outputSubsample, outputSizes[i][1:end-1], λ, concatStart, keeper, j)
          else
            futureThing = @spawn transformMidLayer!(layers, concatOutput, x, dataSizes, outerAxes, innerAxes, innerSub, k, i, daughters, totalScales[i], nonlinear, subsam, outputSubsample, outputSizes[i][1:end-1], concatStart, λ, resultingSize, dataChannel)
          end
        else
          #TODO
          if stType=="decreasing"
            futureThing = @spawn transformFinalLayer!(layers, concatOutput, x, outerAxes, innerAxes, innerSub, k, i, daughters, totalScales[i], nonlinear, subsam, outputSubsample, λ, resultingSize, countInFinalLayer, concatStartLast, keeper, nScalesLastLayer)
          else
            futureThing = @spawn transformFinalLayer!(layers, concatOutput, x, outerAxes, innerAxes, innerSub, k, i, daughters, totalScales[i], nonlinear, subsam, outputSubsample, λ, resultingSize, concatStart, concatStartLast, nScalesLastLayer)
          end
        end
        if typeof(futureThing)<: typeof(Distributed.RemoteException)
          throw(futureThing)
        else
          listOfProcessResults[(j-1)*size(cur[1])[end] + λ] = futureThing
        end
      end
    end

    # using one channel per
    for (λ,x) in enumerate(listOfProcessResults)
      tmpFetch = fetch(x)
      if i<layers.m
        nextData[λ] = tmpFetch
      end
    end
  end
return concatOutput
end

function transformMidLayer!(layers, concatOutput, curPath, dataSize, outerAxes, innerAxes, innerSub, k, i, daughters, nScales, nonlinear, subsam, outputSubsample, outputSizeThisLayer, concatStart, λ, resultingSize, homeChannel)
  # make an array to return the results
  toBeHandedBack = zeros(outputSizeThisLayer..., size(daughters,2)-1)
  # first perform the continuous wavelet transform on the data from the previous layer
  output = cwt(curPath[outerAxes..., innerAxes..., λ], layers.shears[i], daughters, nScales=nScales)
  innerMostAxes = axes(output)[end:end]
  # iterate over the non transformed dimensions of output
  for outer in eachindex(view(output, outerAxes..., [1 for i=1:(k+1)]...))
    # output may have length zero if the starting path has low enough scale
    for j = 2:size(output)[end]
      toBeHandedBack[outer, innerSub..., j-1] = nonlinear.(subsam(output[outer, innerAxes..., j], layers.subsampling[i]))
    end
    # actually write to the output
    tmpOut = Array{Float64,k}(real(subsam(output[outer, innerAxes..., end], layers.subsampling[i])))
    sizeTmpOut = prod(size(tmpOut))
    if outputSubsample[1] > 1 || outputSubsample[2] > 1
      concatOutput[outer, concatStart .+ (λ-1)*resultingSize[i] .+ (0:(resultingSize[i]-1))] = subsam(tmpOut, resultingSize[i], absolute = true)
    else
      concatOutput[outer, concatStart .+ (λ-1)*resultingSize[i] .+ (0:sizeTmpOut-1)] = reshape(tmpOut, (sizeTmpOut))
    end
  end
  return toBeHandedBack
end

function transformFinalLayer!(layers, concatOutput, curPath, outerAxes, innerAxes, innerSub, k, i, daughters, nScales, nonlinear, subsam, outputSubsample, λ, resultingSize, concatStart, concatStartLast, nScalesLastLayer)
  localIndex = 0
  # iterate over the outerAxes
    # first perform the continuous wavelet transform on the data from the previous layer
    output = cwt(curPath[outerAxes..., innerAxes..., λ], layers.shears[i], daughters, nScales=nScales)
  for outer in eachindex(view(output, outerAxes..., [1 for i=1:(k+1)]...))
    # write the output from the m-1st layer to output
    tmpOut = Array{Float64,k}(real(subsam(output[outer, innerAxes..., end], layers.subsampling[i])))
    if outputSubsample[1] > 1 || outputSubsample[2] > 1
      #println("outer, outPos[outer] = $((outer, outPos[outer]))")
      concatOutput[outer, concatStart .+ (λ-1)*resultingSize[i] .+ (0:(resultingSize[i]-1))] = subsam(tmpOut, resultingSize[i], absolute = true)
    else
      sizeTmpOut = prod(size(tmpOut))
      concatOutput[outer, concatStart .+ (λ-1)*resultingSize[i].+(0:(sizeTmpOut-1))] = reshape(tmpOut, (sizeTmpOut))
    end
  end
  # write the output from the mth layer to output
  for j = 2:size(output)[end]
    finalOutput = cwt(curPath[outerAxes..., innerAxes..., λ], layers.shears[i+1], daughters, nScales=1)
    for outer in eachindex(view(finalOutput, outerAxes..., [1 for i=1:(k+1)]...))
      tmpFinal = Array{Float64,k}(real(subsam(finalOutput[outer, innerAxes..., end], layers.subsampling[i])))
      if outputSubsample[1] > 1 || outputSubsample[2] > 1

        concatOutput[outer, concatStartLast + (λ-1)*nScalesLastLayer*resultingSize[i] + (j-2)*resultingSize[i] .+ (0:(resultingSize[i]-1))] = subsam(tmpFinal, resultingSize[i+1], absolute = true)
      else
        sizeTmpOut = prod(size(tmpFinal))
        concatOutput[outer, (concatStartLast + (λ-1)*nScalesLastLayer*resultingSize[i] + localIndex + (λ-1)*resultingSize[i+1]).+(0:(sizeTmpOut-1))] = reshape(tmpFinal, (sizeTmpOut))
      end
    end
  end
end


text:
function hello()
    return "Goodbye, Mars!"
end

ext: jl


text:
"""
    IndexEntry(checkpoint_path, [checkpoint_name, prefixes, tags])

This is an index entry describing the output file from a checkpoint.
You can retrieve a list of these from a folder full of such outputs, using
[`index_checkpoint_files`](@ref).

For accessing details of the IndexEntry the following helpers are provided:
[`checkpoint_path`](@ref), [`checkpoint_name`](@ref), [`prefixes`](@ref), [`tags`](@ref).
Further: `getproperty` is overloaded so that you can access the value of the tag `:foo` via
`x.foo`.
"""
struct IndexEntry
    checkpoint_path::AbstractPath
    checkpoint_name::AbstractString
    prefixes::NTuple{<:Any, AbstractString}
    tags::NTuple{<:Any, Pair{Symbol, <:AbstractString}}
end

IndexEntry(file) = IndexEntry(Path(file))
function IndexEntry(filepath::AbstractPath)
    isfile(filepath) || throw(DomainError(filepath, "Need an existant file."))
    # skip any non-tag directories at the start. Note this will be tricked if those have "="
    # in them but probably not worth handling, unless an issue comes up
    first_tag_ind = something(findfirst(contains("="), filepath.segments), 1)
    segments = filepath.segments[first_tag_ind:end-1]

    prefixes = filter(!contains("="), segments)
    tags = map(filter(contains("="), segments)) do seg
        tag, val = split(seg, "="; limit=2)
        return Symbol(tag)=>val
    end
    checkpoint_name = filename(filepath)
    return IndexEntry(filepath, checkpoint_name, prefixes, tags)
end

"""
    checkpoint_path(x::IndexEntry)

The checkpoint_path to the checkpoint output file.
For example `S3Path("s3::/mybucket/tag1=a/tag2=b/group1/group2/forecasts.jlso)`.
"""
checkpoint_path(x::IndexEntry) = getfield(x, :checkpoint_path)

"""
    checkpoint_name(x::IndexEntry)

The checkpoint_name of the checkpoint output file.
If the checkpoint was saved used `checkpoint(Forecasters, "forecasts", ...)` then it's
`checkpoint_name` is `"forecasts"`.
"""
checkpoint_name(x::IndexEntry) = getfield(x, :checkpoint_name)

"""
    prefixes(x::IndexEntry)

The prefixes of the checkpoint output file. This is a tuple of any prefixes specified.
If the checkpoint was saved used `checkpoint(Forecasters, "forecasts", ...)` then it's
`prefixes` are `("Forecasters",)`. Generally in practice there are either one matching the
module the checkpoint was declared in (as in the previous example), or zero if it wasn't
specified (as in `checkpoint("forecasts",...)`). If it was saved as
`checkpoint("Foo.Bar", "forecasts.jlso")` then prefixed will return `("Foo", "Bar")`.
"""
prefixes(x::IndexEntry) = getfield(x, :prefixes)

"""
    tags(x::IndexEntry)

All tags and there values of the checkpoint. This is a collection of pairs.
For example if the checkpoint was saved using:
```julia
with_tag(:sim_now=>DateTime(2000,1,1,9,00), grid=ERCOT) do
    checkpoint(Forecasters, "forecasts", ...)`
end
```
Then `tags(x)` would return: `(:sim_now="2000-01-01T09:00:00", :grid=>"ERCOT")`.

Note that if the tags are unique, then their values call also be accessed via a
`getproperty` overload, e.g as `x.sim_now` or `x.grid`.
"""
tags(x::IndexEntry) = getfield(x, :tags)

_tag_names(x::IndexEntry) = first.(tags(x))

#Tables.columnnames(x::IndexEntry) = propertynames(x)
function Base.propertynames(x::IndexEntry)
    return [:prefixes, :checkpoint_name, _tag_names(x)..., :checkpoint_path]
end

function Base.getproperty(x::IndexEntry, name::Symbol)
    inds = findall(==(name)∘first, tags(x))
    if length(inds) > 1
        error(
            "The checkpoint was tagged with $name multiple times (positions $inds). " *
            "you thus can not use the getproperty shorthand to check its tag, as we don't" *
            "know which you mean. Use the `tags` function instead."
        )
    elseif length(inds) == 0
        if hasfield(IndexEntry, name)
            return getfield(x, name)
        else
            error(
                "The checkpoint does not have the tag $name. It has tags: $(_tag_names(x))"
            )
        end
    else  # tag found
        if hasfield(IndexEntry, name)  #also is a field, so lets error
            error(
                "$name is both a tag, and a field of the index entry itself. use the " *
                "function $name(x) to get the index entry field, or work with `tags(x)` " *
                "for the field"
            )
        else
            _, val = only(tags(x)[inds])
            return val
        end
    end
end

"""
    index_checkpoint_files(dir)

Constructs a index for all the files output by checkpoints located within  `dir`.
This index tells you their checkpoint_name, checkpoint_path, tags, etc.
See [`IndexEntry`](@ref) for full information on what is recorded.

Handily, a `Vector{IndexEntry}` is a valid [Tables.jl Table][1]. This means you can
do `DataFrame(index_checkpoint_files(dir))` and get a nice easy to work with [DataFrame][2].

You can also work with it directly, say you wanted to get all checkpoints files for
`forecasts` with the tag `model="Stage1"`. This would be something like:
```julia
[checkpoint_path(x) for x in index_checkpoint_files(dir) if x.checkpoint_name=="forecast" && x.model=="Stage1"]
```

1: https://github.com/JuliaData/Tables.jl
2: https://github.com/JuliaData/DataFrames.jl
"""
function index_checkpoint_files(dir::AbstractPath)
    map(Iterators.filter(==("jlso") ∘ extension, walkpath(dir))) do checkpoint_path
        return IndexEntry(checkpoint_path)
    end
end

index_checkpoint_files(dir) = index_checkpoint_files(Path(dir))

"""
    index_files(dir)

Constructs a index for all the files located within  `dir`.
Same as [`index_checkpoint_files`] except not restricted to files created by Checkpoints.jl.
"""
function index_files(dir::AbstractPath)
    map(Iterators.filter(isfile, walkpath(dir))) do path
        return IndexEntry(path)
    end
end

index_files(dir) = index_files(Path(dir))


text:
module Model_Household

  export Household, Person, People_Dict
  export uprate!, equivalence_scale, oldest_person, default_bu_allocation
  export get_benefit_units, num_people, get_head,get_spouse, printpids

  include( "household.jl" )

end


text:
# cross-correlation module
export clean_up!, clean_up, correlate, compute_cc, generate_pairs, map_cc, corrmap
export coherence!, coherence, deconvolution!, deconvolution, whiten, whiten!

"""
    clean_up!(A,freqmin,freqmax,fs)

Demean, detrend, taper and filter time series.

# Arguments
- `A::AbstractArray`: Time series.
- `fs::Float64`: Sampling rate of time series `A` in Hz.
- `freqmin::Float64`: Pass band low corner frequency in Hz.
- `freqmax::Float64`: Pass band high corner frequency in Hz.
"""
function clean_up!(A::AbstractArray, freqmin::Real, freqmax::Real, fs::Real;
                   corners::Int=4, zerophase::Bool=true)
    demean!(A)
    detrend!(A)
    taper!(A,fs)
    bandpass!(A,freqmin,freqmax,fs,corners=corners,zerophase=zerophase)
    return nothing
end
clean_up(A::AbstractArray, freqmin::Real, freqmax::Real, fs::Real;
         corners::Int=4, zerophase::Bool=false) =
                 (U = deepcopy(A); clean_up!(U,freqmin,freqmax, fs,
                  corners=corners, zerophase=zerophase); return U)

clean_up!(C::CorrData,freqmin::Float64,freqmax::Float64; corners::Int=4,
          zerophase::Bool=true) = clean_up!(C.corr,freqmin,freqmax,C.fs,
          corners=corners,zerophase=zerophase)

clean_up!(R::RawData,freqmin::Float64,freqmax::Float64; corners::Int=4,
          zerophase::Bool=true) = clean_up!(C.x,freqmin,freqmax,R.fs,
          corners=corners,zerophase=true)

"""
    correlate(FFT1, FFT2, N, maxlag, corr_type='cross-correlation')

Cross-correlate ambient noise data in the frequency domain.

Cross-correlation can be done using one of three options:

- Cross-correlation: ``C_{AB}(ω) = u_A(ω) u^∗_B(ω)``
- Coherence: ``C_{AB}(ω) = \frac{u_A(ω) u^∗_B(ω)}{∣ u_A(ω) ∣  ∣ u_B(ω) ∣}``
- Deconvolution: ``C_{AB}(ω) = \frac{u_A(ω) u^∗_B(ω)}{∣ u_B(ω) ∣^2}``

Smoothing of FFTs for coherence and deconvolution should be done before
cross-correlating.

# Arguments
- `FFT1::AbstractArray`: Complex Array of fourier transform of ambient noise data.
- `FFT2::AbstractArray`: Complex Array of fourier transform of ambient noise data.
- `N::Int`: Number of input data points in time domain, equal to `cc_len` * `fs`.
- `maxlag::Int`: Number of data points in cross-correlation to save,
                 e.g. `maxlag = 2000` will save lag times = -2000/fs:2000/fs s.
- `corr_type::String`: Type of correlation: `cross-correlation`, `coherence` or
                       `deconv`.
"""
function correlate(FFT1::AbstractArray, FFT2::AbstractArray, N::Int, maxlag::Int)

    corrF = conj.(FFT1) .* FFT2

    # take inverse fft
    corrT = irfft(corrF,N,1)
    corrT = fftshift(corrT,1)

    # return corr[-maxlag:maxlag]
    t = range(-Int(N/2) + 1, stop=Int(N/2) - 1)
    ind = findall(abs.(t) .<= maxlag)
    corrT = corrT[ind,:]
end

"""
    compute_cc(FFT1::FFTData, FFT2::FFTData, maxlag::Float64;
               corr_type::String="cross-correlation")

Cross-correlate ambient noise data in the frequency domain.

Cross-correlation can be done using one of three options:

- Cross-correlation: ``C_{AB}(ω) = u_A(ω) u^∗_B(ω)``
- Coherence: ``C_{AB}(ω) = \frac{u_A(ω) u^∗_B(ω)}{∣ u_A(ω) ∣  ∣ u_B(ω) ∣}``
- Deconvolution: ``C_{AB}(ω) = \frac{u_A(ω) u^∗_B(ω)}{∣ u_B(ω) ∣^2}``

# Arguments
- `FFT1::FFTData`: FFTData object of fft'd ambient noise data.
- `FFT2::FFTData`: FFTData object of fft'd ambient noise data.
- `maxlag::Float64`: Maximum lag time (in seconds) in cross-correlation to save,
                     e.g. `maxlag = 20.` will save lag times = -20.:20. s.
- `corr_type::String`: Type of correlation: `cross-correlation`, `coherence` or
                       `deconv`.
"""
function compute_cc(FFT1::FFTData, FFT2::FFTData, maxlag::Float64;
                    corr_type::String="cross-correlation")

    N = convert(Int,round(FFT1.cc_len * FFT1.fs)) # number of data points
    comp = FFT1.name[end] * FFT2.name[end]
    # get intersect of dates; return nothing if no intersect
    inter = intersect(FFT1.t,FFT2.t)
    if length(inter) == 0
        return nothing
    end

    ind1 = findall(x -> x ∈ inter, FFT1.t)
    ind2 = findall(x -> x ∈ inter, FFT2.t)

    corr = correlate(FFT1.fft[:,ind1], FFT2.fft[:,ind2], N,
                     convert(Int,round(maxlag * FFT1.fs)))
    rotated = false

    return CorrData(FFT1, FFT2, comp, rotated, corr_type,
                    maxlag, inter, corr)

end

"""
    corrmap(A,maxlag,smoothing_half_win,corr_type,OUTDIR)

Compute cross-correlations using a parallel map.

`corrmap` takes in an array of FFT's and cross-correlates using full parallism.
For a list of `N` `FFT`s, there are `N * (N - 1) / 2` total possible
cross-correlations. To avoid parallel I/O issues, `corrmap` copies the current
correlation `F`, then maps the correlation of `F` against all other correlations
in parallel.

# Arguments
- `A::AbstractArray`: Array of FFTData objects.
- `maxlag::Float64`: Maximum lag time (in seconds) in cross-correlation to save,
                     e.g. `maxlag = 20.` will save lag times = -20.:20. s.
- `smoothing_half_win::Int`: Number of points to smooth spectrum of
                    cross-correlations if using deconvolution or coherence.
- `corr_type::String`: Type of correlation: `cross-correlation`, `coherence` or
                       `deconv`.
- `OUTDIR::String`: Path to save correlation, e.g. "/home/ubuntu/CORR/".
"""
function corrmap(A::Array{FFTData,1},maxlag::Float64, smoothing_half_win::Int,
                 corr_type::String,OUTDIR::String;
                 water_level::Union{Nothing,AbstractFloat}=nothing)
        N = size(A,1)

        if corr_type == "coherence"
            A = pmap(coherence,A,fill(smoothing_half_win,N),fill(water_level,N))
        end

        # copy the current FFT and correlate against all remaining
        for ii = 1:N-1
            FFT = pop!(A)
            if corr_type == "deconv"
                deconvolution!(FFT,smoothing_half_win,water_level)
            end
            out = fill(deepcopy(FFT),length(A))
            pmap(map_cc,out,A,fill(maxlag,length(A)),
                 fill(corr_type,length(A)),fill(OUTDIR,length(A)))
        end
end

function corrmap(A::Array{FFTData,1},maxlag::Float64, smoothing_half_win::Int,
                 corr_type::String,OUTDIR::String,
                 interval::Union{Month,Day,Hour,Second};
                 water_level::Union{Nothing,AbstractFloat}=nothing)
        N = size(A,1)

        if corr_type == "coherence"
            A = pmap(coherence,A,fill(smoothing_half_win,N),fill(water_level,N))
        end

        # copy the current FFT and correlate against all remaining
        for ii = 1:N-1
            FFT = pop!(A)
            out = fill(deepcopy(FFT),length(A))
            if corr_type == "deconv"
                deconvolution!(FFT,smoothing_half_win,water_level)
            end
            pmap(map_cc,out,A,fill(maxlag,length(A)),
                 fill(corr_type,length(A)),fill(OUTDIR,length(A)),
                 fill(interval,length(A)))
        end
end

"""
  map_cc(FFT1,FFT2,maxlag,OUTDIR)

Input function for corrmap.

Correlates `FFT1` and `FFT2`.

# Arguments
- `FFT1::FFTData`: FFTData object of fft'd ambient noise data.
- `FFT2::FFTData`: FFTData object of fft'd ambient noise data.
- `maxlag::Float64`: Maximum lag time (in seconds) in cross-correlation to save,
                 e.g. `maxlag = 20.` will save lag times = -20.:20. s.
- `corr_type::String`: Type of correlation: `cross-correlation`, `coherence` or
                   `deconv`.
- `OUTDIR::String`: Path to save correlation, e.g. "/home/ubuntu/CORR/".

"""
function map_cc(FFT1::FFTData,FFT2::FFTData,maxlag::Float64,
                corr_type::String,OUTDIR::String)
    println("Correlation $(FFT1.name), $(FFT2.name)")
    C = compute_cc(FFT1,FFT2,maxlag,corr_type=corr_type)
    if !isnothing(C)
        save_corr(C,OUTDIR)
    end
    return nothing
end

function map_cc(FFT1::FFTData,FFT2::FFTData,maxlag::Float64,
                corr_type::String,OUTDIR::String,
                interval::Union{Month,Day,Hour,Second})
    println("Correlation $(FFT1.name), $(FFT2.name)")
    C = compute_cc(FFT1,FFT2,maxlag,corr_type=corr_type)
    if !isnothing(C)
        stack!(C,interval=interval)
        save_corr(C,OUTDIR)
    end
    return nothing
end

"""
    whiten(A, freqmin, freqmax, fs, pad=50)

Whiten spectrum of rfft `A` between frequencies `freqmin` and `freqmax`.
Returns the whitened rfft of the time series.

# Arguments
- `A::AbstractArray`: Time series.
- `fs::Real`: Sampling rate of time series `A`.
- `freqmin::Real`: Pass band low corner frequency.
- `freqmax::Real`: Pass band high corner frequency.
- `N::Int`: Number of input time domain samples for each rfft.
- `pad::Int`: Number of tapering points outside whitening band.
"""
function whiten!(A::AbstractArray, freqmin::Float64, freqmax::Float64, fs::Float64,
                 N::Int;pad::Int=50)
    T = real(eltype(A))

    # get size and convert to Float32
    Nrows,Ncols = size(A)

    # get whitening frequencies
    freqvec = FFTW.rfftfreq(N,fs)
    left = findfirst(x -> x >= freqmin, freqvec)
    right = findlast(x -> x <= freqmax, freqvec)
    low, high = left - pad, right + pad

    if low <= 1
        low = 1
        left = low + pad
    end

    if high > length(freqvec)
        high = length(freqvec)- 1
        right = high - pad
    end

    # left zero cut-off
    compzero = complex(T(0))
    for jj = 1:Ncols
         for ii = 1:low
            A[ii,jj] = compzero
        end

    # left tapering
         for ii = low+1:left
            A[ii,jj] = cos(T(pi) ./ T(2) .+ T(pi) ./ T(2) .* (ii-low-1) ./ pad).^2 .* exp(
            im .* angle(A[ii,jj]))
        end

    # pass band
         for ii = left+1:right
            A[ii,jj] = exp(im .* angle(A[ii,jj]))
        end

    # right tapering
         for ii = right+1:high
            A[ii,jj] = cos(T(pi) ./ T(2) .* (ii-right) ./ pad).^2 .* exp(
            im .* angle(A[ii,jj]))
        end

    # right zero cut-off
         for ii = high+1:size(A,1)
            A[ii,jj] = compzero
        end
    end
    return nothing
end
whiten(A::AbstractArray, freqmin::Float64, freqmax::Float64, fs::Float64, N::Int;
       pad::Int=50) = (U = deepcopy(A);
                       whiten!(U,freqmin,freqmax,fs,N,pad=pad);
                       return U)

"""
   whiten(F, freqmin, freqmax)

Whiten spectrum of FFTData `F` between frequencies `freqmin` and `freqmax`.
Uses real fft to speed up computation.
Returns the whitened (single-sided) fft of the time series.

# Arguments
- `F::FFTData`: FFTData object of fft'd ambient noise data.
- `freqmin::Real`: Pass band low corner frequency.
- `freqmax::Real`: Pass band high corner frequency.
- `pad::Int`: Number of tapering points outside whitening band.
"""
function whiten!(F::FFTData, freqmin::Float64, freqmax::Float64;pad::Int=50)
    if freqmin < F.freqmin && freqmax > F.freqmax
        @warn "Whitening frequencies ($freqmin, $freqmax Hz) are wider than frequencies
        in FFTData ($(F.freqmin),$(F.freqmax) Hz). Whitening in ($(F.freqmin),$(F.freqmax) Hz) band."
    elseif freqmin < F.freqmin
        @warn "Low whitening frequency $freqmin Hz is lower than minumum frequency
        in FFTData ($(F.freqmin) Hz). Whitening in ($(F.freqmin),$freqmax Hz) band."
    elseif freqmax > F.freqmax
        @warn "High whitening frequency $freqmax Hz is higher than maximum frequency
        in FFTData ($(F.freqmax) Hz). Whitening in ($freqmin,$(F.freqmax) Hz) band."
    end

    freqmin = max(freqmin,F.freqmin) # check for freqmin = 0
    freqmax = min(freqmax,max(F.freqmax,1 / F.cc_len)) # check for freqmax = 0
    whiten!(F.fft, freqmin, freqmax, F.fs,convert(Int, F.fs * F.cc_len), pad=pad)
    F.whitened = true
    F.freqmin = freqmin
    F.freqmax = freqmax
    return nothing
end
whiten(F::FFTData, freqmin::Float64, freqmax::Float64;pad::Int=50) =
      (U = deepcopy(F); whiten!(U,freqmin,freqmax,pad=pad);return U)

function whiten!(R::RawData,freqmin::Float64, freqmax::Float64; pad::Int=50)
    if freqmin < R.freqmin && freqmax > R.freqmax
        @warn "Whitening frequencies ($freqmin, $freqmax Hz) are wider than frequencies
        in RawData ($(R.freqmin),$(R.freqmax) Hz). Whitening in ($(R.freqmin),$(R.freqmax) Hz) band."
    elseif freqmin < R.freqmin
        @warn "Low whitening frequency $freqmin Hz is lower than minumum frequency
        in RawData ($(R.freqmin) Hz). Whitening in ($(R.freqmin),$freqmax Hz) band."
    elseif freqmax > R.freqmax
        @warn "High whitening frequency $freqmax Hz is higher than maximum frequency
        in FFTData ($(R.freqmax) Hz). Whitening in ($freqmin,$(R.freqmax) Hz) band."
    end

    freqmin = max(freqmin,R.freqmin) # check for freqmin = 0
    freqmax = min(freqmax,max(R.freqmax,1 / R.cc_len)) # check for freqmax = 0
    FFT = rfft(R.x,1)
    whiten!(FFT,freqmin,freqmax,R.fs, convert(Int, R.fs * R.cc_len), pad=pad)
    R.x .= irfft(FFT,convert(Int,R.cc_len * R.fs),1)
    R.freqmin = freqmin
    R.freqmax = freqmax
    R.whitened = true
    return nothing
end
whiten(R::RawData,freqmin::Float64, freqmax::Float64; pad::Int=50) =
      (U = deepcopy(R); whiten!(U,freqmin,freqmax,pad=pad);return U)

"""

  coherence!(F,half_win, water_level)

Apply coherence method to FFTData `F`. Where,
``C_{AB}(ω) = \frac{u_A(ω) u^∗_B(ω)}{∣ u_A(ω) ∣  ∣ u_B(ω) ∣}``

# Arguments
- `F::FFTData`: FFTData object of fft'd ambient noise data.
- `half_win::Int`: Number of points in half-window to smooth spectrum.
- `water_level::AbstractFloat`: Regularization parameter for spectral smoothing.
                                0.01 is a common value [Mehta, 2007].
"""
function coherence!(F::FFTData, half_win::Int,
                    water_level::Union{Nothing,AbstractFloat}=nothing)
    smoothF = smooth(abs.(F.fft),half_win)
    if !isnothing(water_level)
        reg = water_level .* mean(abs.(F.fft),dims=1)
        smoothF .+= reg
    end
    F.fft ./= smoothF
end
coherence(F::FFTData,half_win::Int,
          water_level::Union{Nothing,AbstractFloat}=nothing) =
          (U = deepcopy(F);coherence!(U,half_win,water_level);
          return U)

"""

  deconvolution!(F,half_win, water_level)

Apply deconvolution method to FFTData `F`. Where,
``C_{AB}(ω) = \frac{u_A(ω) u^∗_B(ω)}{∣ u_B(ω) ∣^2}``

# Arguments
- `F::FFTData`: FFTData object of fft'd ambient noise data.
- `half_win::Int`: Number of points in half-window to smooth spectrum.
- `water_level::AbstractFloat`: Regularization parameter for spectral smoothing.
                              0.01 is a common value [Mehta, 2007].
"""
function deconvolution!(F::FFTData, half_win::Int,
                        water_level::Union{Nothing,AbstractFloat}=nothing)
    smoothF = smooth(abs.(F.fft).^2,half_win)
    if !isnothing(water_level)
        reg = water_level .* mean(abs.(F.fft).^2,dims=1)
        smoothF .+= reg
    end
    F.fft ./= smoothF
end
deconvolution(F::FFTData,half_win::Int,
              water_level::Union{Nothing,AbstractFloat}=nothing) =
              (U = deepcopy(F);deconvolution!(U,half_win,water_level);
              return U)

function generate_pairs(files::AbstractArray)
    N = length(files)
    num_pairs = convert(Int,round(N * (N-1) / 2 + N))
    pairs = Array{Array{String,1},1}(undef,num_pairs)
    count = 0
    for ii = 1:length(files)
      for jj = ii:length(files)
          count += 1
          pairs[count] = [files[ii], files[jj]]
      end
    end
    return pairs
end


text:
# Run advection tests

include(joinpath(Pkg.dir("PDESolver"), "src/input/make_input.jl"))

using PDESolver
#using Base.Test
using FactCheck
using ODLCommonTools
using PdePumiInterface  # common mesh interface - pumi
using SummationByParts  # SBP operators
using AdvectionEquationMod
using NonlinearSolvers   # non-linear solvers
using ArrayViews

fname = "input_vals_newton.jl"

@testset "----- Testing Newtons Method: Finite Difference -----" begin
  mesh, sbp, eqn, opts = solvePDE(fname)
  @test  calcNorm(eqn, eqn.res_vec, strongres=true)  < opts["res_abstol"]
end

@testset "----- Testing Newtons Method: Complex Step -----" begin
  arg_dict["jac_method"] = 2
  make_input(arg_dict, "input_vals_newton2")
  fname = "input_vals_newton2.jl"
  mesh, sbp, eqn, opts = solvePDE(fname)
  @test  calcNorm(eqn, eqn.res_vec, strongres=true)  < opts["res_abstol"]
end

@testset "----- Testing Newtons Method: Petsc -----" begin
  arg_dict["jac_type"] = 3
  make_input(arg_dict, "input_vals_newton3")
  fname = "input_vals_newton3.jl"
  mesh, sbp, eqn, opts = solvePDE(fname)
  @test  calcNorm(eqn, eqn.res_vec, strongres=true)  < opts["res_abstol"]
end


text:
@testset "WeightedSelector" begin
    Random.seed!(123)
    s = WeightedSelector(true)
    values = [0.2, 0.3, 0.5]
    action_counts = countmap([s(values; step = i) for i = 1:1000])
    for i = 1:length(values)
        @test isapprox(action_counts[i] / 1000, values[i], atol = 0.1)
    end
end


text:
function b4step1(mbc::Int, mx::Int, meqn::Int, q, xloweri::Float64, dxi::Float64, t::Float64, dt::Float64, mauxi::Int, aux)
#
#    ! Called before each call to step1.
#    ! Use to set time-dependent aux arrays or perform other tasks.
#    !
#    ! This default version does nothing.
#
#    implicit none
#    integer, intent(in) :: mbc,mx,meqn,maux
#    real(kind=8), intent(in) :: xlower,dx,t,dt
#    real(kind=8), intent(inout) :: q(meqn,1-mbc:mx+mbc)
#    real(kind=8), intent(inout) :: aux(maux,1-mbc:mx+mbc)
#
end # subroutine b4step1


text:
# default options go here
# if user specifies option they will be overwritten
export CurrentOptions
export InvalidFermiOption
export @reset
export @set
export @get
export @molecule

const DefaultOptions = Dict{String,Union{Float64,Int,String,Bool,Nothing}}(
                                  "molstring" => """
                                  O        1.2091536548      1.7664118189     -0.0171613972
                                  H        2.1984800075      1.7977100627      0.0121161719
                                  H        0.9197881882      2.4580185570      0.6297938832
                                  """,
                                  "basis" => "sto-3g",
                                  "jkfit" => "auto",
                                  "rifit" => "auto",
                                  "charge" => 0,
                                  "multiplicity" => 1,
                                  "unit" => "angstrom",
                                  "reference" => "rhf",
                                  "scf_max_iter" => 50,
                                  "scf_max_rms" => 10^-10,
                                  "scf_alg" => "conventional",
                                  "oda" => true,
                                  "oda_cutoff" => 1E-1,
                                  "oda_shutoff" => 20,
                                  "scf_guess" => "gwh",
                                  "quiet" => true,
                                  "mp2_type" => "df",
                                  "precision" => "double",
                                  "cc_alg" => "CTF",
                                  "ci_alg" => "aci",
                                  "σ"      => 0.001,
                                  "γ"      => 1.0,
                                  "e_conv" => 10,
                                  "d_conv" => 8,
                                  "cc_max_iter" => 50,
                                  "cc_max_rms" => 10^-10,
                                  "cc_e_conv" => 10^-10,
                                  "bcc_max_t1" => 1^-7,
                                  "preconv_T1" => true,
                                  "drop_occ" => 0,
                                  "drop_vir" => 0,
                                  "diis" => true,
                                  "cc_diis" => true,
                                  "ndiis" => 8,
                                  "cc_ndiis" => 3,
                                  "bcc_tol" => 1.0e-6,
                                  "diis_prec" => "single",
                                  "diis_start" => 3,
                                  "cc_damp_ratio" => 0.0,
                                  "cc_diis_relax" => 3,
                                  "num_frozen" => 0,
                                  "aci_print_screen" => nothing,
                                  "cas_frozen" => 0,
                                  "cas_active" => -1,
                                  "cas_cutoff" => 10^-9,
                                  "cas_nroot" => 1,
                                  "min_matrix_elem" => 10^-9,
                                  "precision_override" => false
                                 )
"""
    Fermi.CurrentOptions

Dictionary containing options for Fermi. Any information not given
explicitly to Methods is obtained from here.
"""
CurrentOptions = deepcopy(DefaultOptions)

"""
    Fermi.@reset(key=nothing)

Restores `key` option to default setting. If `key == nothing`, restores all options to their default values.
"""
macro reset(key=nothing)
    if key == nothing
        return quote
            for key in keys(Fermi.CurrentOptions)
                Fermi.CurrentOptions[key] = Fermi.DefaultOptions[key]
            end
        end
    elseif String(key) in keys(Fermi.CurrentOptions)
        return quote
            Fermi.CurrentOptions[$(String(key))] = Fermi.DefaultOptions[$(String(key))]
        end
    else
        return quote
            error("Invalid option "*$(String(key)))
        end
    end
end

"""
    Fermi.@set

Set options for Fermi computations. It saves the options into Fermi.CurrentOptions.

*Usage:*  @set A B

A is set to B. By default A is taken as a string. B is evaluated at runtime. If the evaluation
is not possible, B is converted to a string.

# Examples

```
@set basis cc-pVDZ
@set cc_max_iter 100
@set e_conv 10^-9
```
One can also use the block syntax
```
@set {
    basis cc-pVDZ
    cc_max_iter 100
    e_conv 10^-9
}
```
Note that for the block syntax variables are not accepted because the evaluation of B is
done at parse time.
```
mybasis = 6-31g
@set {
    basis mybasis
}
```
Will set the basis to "mybasis" not "6-31g".
"""
macro set(opt,val)
    clean_up(s) = String(filter(c->!occursin(c," ():"),s))
    A = clean_up(repr(opt))
    B = clean_up(repr(val))
    quote
        try
            CurrentOptions[$A] = $val
        catch UndefVarError
            CurrentOptions[$A] = $B
        end
    end |> esc
end

macro set(block)
    lines = split(repr(block),";")
    quote
        clean_up(s) = String(strip(filter(c->!occursin(c," {}():"),s)))
        for l in $lines
            opt, val = split(strip(l), " ", limit=2)
            opt = clean_up(opt)
            val = clean_up(val)
            try
                CurrentOptions[opt] = eval(Meta.parse(val))
            catch UndefVarError
                CurrentOptions[opt] = val
            end
        end
    end
end

macro get(opt)
    clean_up(s) = String(filter(c->!occursin(c," ():"),s))
    A = clean_up(repr(opt))
    quote
        try
            CurrentOptions[$A]
        catch KeyError
            throw(InvalidFermiOption($A*" not a valid option."))
        end
    end |> esc
end
"""
    Fermi.@molecule

Set the molecule used in computations. String is save into Fermi.CurrentOptions

# Example
```
@molecule {
C                 -0.00000000     0.00000000    -0.00000000
O                  0.00000000     0.00000000    -2.19732722
O                 -0.00000000    -0.00000000     2.19732722
}
```
"""
macro molecule(block)
    clean_up(s) = strip(filter(c->!occursin(c,"{}():"),s))
    mol = repr(block)
    mol = replace(mol, ";"=>"\n")
    mol = clean_up(mol)
    CurrentOptions["molstring"] = String(mol)
end

struct InvalidFermiOption <: Exception
    msg::String
end
Base.showerror(io::IO, e::InvalidFermiOption) = print(io, "InvalidFermiOption: ", e.msg)

#function notimplemented()
#    @output "🚧 Not implemented yet! We're working on it 🔨 👷 \n"
#end


text:
import EnergyExpressions: BitPattern

@testset "Show binary values" begin
    @test string(BitPattern(1)) == "1"
    @test string(BitPattern(2)) == ".1"
    @test string(BitPattern(3)) == "11"
    @test string(BitPattern(5)) == "1.1"
    @test string(BitPattern(16)) == ".... 1"
end


text:
@testset "Conversions" begin
    x = 1.25
    i = 4

    # Floats test conversion back and forth
    @test x == Float64(Posit8(x))
    @test x == Float64(Posit16(x))
    @test x == Float64(Posit32(x))

    @test x == Float32(Posit8_2(x))
    @test x == Float32(Posit16_2(x))
    @test x == Float32(Posit24_2(x))

    # Posit to Posit back and forth
    if Sys.isapple() || Sys.iswindows()
        @test_skip Posit8(x) == Posit8(Posit16(x))
        @test_skip Posit32(x) == Posit32(Posit16(Posit8_2(x)))
    else
        @test Posit8(x) == Posit8(Posit16(x))
        @test Posit32(x) == Posit32(Posit16(Posit8_2(x)))
    end

    # Int to posit back and forth
    @test i == Int64(Posit8(i))
    @test i == Int64(Posit16(i))
    @test i == Int64(Posit32(i))

    # Boolean tests
    @testset for T in (Posit8,Posit16,Posit32)
        T(x) == true*T(x)
        zero(T) == false*T(x)
        T(x)+one(T) == true + T(x)
        T(x) == false + T(x)
        T(true) == one(T)
        T(false) == zero(T)
    end
end


text:
module aoc_24
    export read_input, part1, part2, test_part1, test_part2, factorial
    export direction_dict

    using OffsetArrays, ImageFiltering, Test

    const directions_r = r"e|se|sw|w|nw|ne"

    "Map hexagonal neighbors to cartesian grid"
    const direction_dict = Dict("e" => (0,1), "se" => (1,0), "sw" => (1,-1), "w" => (0,-1), "nw" => (-1,0), "ne" => (-1,1))

    function read_input(filename = "input.txt")
        out = [[m.match for m in eachmatch(directions_r, line)] for line in eachline(filename)]
    end

    const input = read_input()
    const demo = read_input("demo.txt")
    const mini_demo = read_input("demo2.txt")

    function part1(input = input)
        A = parse(input)
        sum(A)
        #282
    end

    "Creates a cartesian image grid"
    function parse(directions = input)
        offsets = [mapreduce(x->aoc_24.direction_dict[x], (x,y)->x .+ y, d) for d in directions]
        ci = CartesianIndex.(offsets)
        ci_max = maximum(ci)
        ci_min = minimum(ci)
        sz = ci_max - ci_min + CartesianIndex(1,1)
        A = OffsetArray( falses(sz.I), OffsetArrays.Origin(ci_min.I...) )
        for c in ci
            A[c] = ~A[c]
        end
        A
    end

    "For Part 2"
    function run( input = input , steps = 100)
        F = OffsetArray(zeros(3,3), -1:1, -1:1)
        F[CartesianIndex.(values(direction_dict))] .= 1
        F
        A = parse(input)

        max_steps = maximum( steps )

        black_tiles = Vector{Int}()

        for step=1:max_steps
            A = padarray(A, Fill(0,(1,1)))

            neighbors = imfilter(A, F, Fill(0))

            # Any black tile with zero or more than 2 black tiles immediately adjacent to it is flipped to white.
            # Any white tile with exactly 2 black tiles immediately adjacent to it is flipped to black.
            flip = ( A .& ( (neighbors .> 2) .| (neighbors .== 0) ) ) .| ( .~A .& (neighbors .== 2) )
            A[ flip ] = .~A[ flip ]

            A = compact(A)

            if step in steps
                push!(black_tiles, sum(A) )
            end

        end
        black_tiles, A
    end

    function part2()
        run()[1][1]
        # 3445
    end

    function test_part1()
        @testset "Part 1" begin
            @test part1(mini_demo) == 2
            @test part1(demo) == 10
            @test part1() == 282
        end
    end
    function test_part2()
        @testset "Part 2" begin
            @test run( demo, [ 1:10 ; 20:10:100 ] )[1] ==[15, 12, 25, 14, 23, 28, 41, 37, 49, 37, 132, 259, 406, 566, 788, 1106, 1373, 1844, 2208]
            @test part2() == 3445
        end
    end

    function test()
        test_part1()
        test_part2()
    end

    "Shrink cartesian image to only contain area that has black tiles"
    function compact(universe)
        # Originall from aoc_17
        e = extrema(findall(universe))
        v = view(universe, [a:b for (a,b) in zip(Tuple.(e)...)]...)
        copy(v)
    end

end


text:
"""
    strength(plate[, s])

Compute the bound vortex sheet strength of the plate.

`s` can be either a single normalized arc length coordinate (between
-1 and 1), or a whole array of coordinates.
"""
strength(plate) = strength(plate, plate.ss)

function strength(plate::Plate{T}, s::Real) where {T}
    @get plate (A, α, B₀, B₁, L, Γ, N)

    γ = zero(T)

    U₋ = 1.0
    U  = 2s
    for n in 2:N-1
        γ += 2A[n]*U
        U₋, U = U, 2s*U - U₋
    end
    γ *= (s^2 - 1)

    γ += 2Γ/(L*π)
    γ += 2(A[0] - B₀)*s
    γ +=  (A[1] - B₁)*(2s^2 - 1)

    if abs2(γ) < 256eps()
        γ = 0.0
    else
        γ /= √(1 - s^2)
    end
    return γ
end

function strength(plate::Plate{R}, ss::AbstractArray{T}) where {R, T <: Real}
    γs = zeros(R, length(ss))
    strength!(γs, plate, ss)
    return γs
end

"""
    strength!(γs, plate[, ss])

Compute the bound vortex sheet strength of `plate` and store it in `γs`.

If an array, `ss`, with normalized arc length coordinates is omitted,
then the strength will be computed at the plate's Chebyshev nodes.
"""
function strength!(γs, plate, ss = plate.ss)
    for (i, s) in enumerate(ss)
        γs[i] = strength(plate, s)
    end
    nothing
end

"""
    bound_circulation(plate[, s])

Compute the bound circulation between the trailing edge of the plate to `s`.

`s` can be either a single normalized arc length coordinate (between
-1 and 1), or a whole array of coordinates.
"""
bound_circulation(plate) = bound_circulation(plate, plate.ss)

function bound_circulation(plate, s::Real)
    @get plate (A, B₀, B₁, L, Γ)
    _bound_circulation(A, B₀, B₁, L, Γ, s)
end

function _bound_circulation(A, B₀, B₁, L, Γ, s)
    N = length(A)

    Γₛ = zero(eltype(A))

    Γₛ  = 2(A[0] - B₀)
    Γₛ +=  (A[1] - B₁)*s

    U₋₂ = 1.0
    U₋  = 2s
    U   = 4s^2 - 1

    for n in 2:N-1
        Γₛ += A[n]*(U/(n+1) - U₋₂/(n-1))
        U₋₂, U₋, U = U₋, U, 2s*U - U₋
    end
    Γₛ *= -0.5L*√(1 - s^2)

    Γₛ -= Γ*(acos(s)/π - 1)

    return Γₛ
end

function bound_circulation(plate::Plate{R}, ss::AbstractArray{T}) where {R,T <: Real}
    Γs = zeros(R, length(ss))
    bound_circulation!(Γs, plate, ss)
    return Γs
end

"""
    bound_circulation!(Γs, plate[, ss])

Compute the bound circulation between the trailing edge of the plate to `ss`, then store it in `Γs`.

If an array, `ss`, with normalized arc length coordinates is omitted,
then the circulation will be computed at the plate's Chebyshev nodes.
"""
function bound_circulation!(Γs, plate, ss = plate.ss)
    for (i, s) in enumerate(ss)
        Γs[i] = bound_circulation(plate, s)
    end
    nothing
end


text:
# === Below Generated by make_jni2.jl ===

get_version(penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetVersion, jint, (Ptr{JNIEnv},), penv)

define_class(name::AbstractString, loader::jobject, buf::Array{jbyte,1}, len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].DefineClass, jclass, (Ptr{JNIEnv}, Cstring, jobject, Ptr{jbyte}, jsize,), penv, name, loader, buf, len)

find_class(name::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].FindClass, jclass, (Ptr{JNIEnv}, Cstring,), penv, name)

from_reflected_method(method::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].FromReflectedMethod, jmethodID, (Ptr{JNIEnv}, jobject,), penv, method)

from_reflected_field(field::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].FromReflectedField, jfieldID, (Ptr{JNIEnv}, jobject,), penv, field)

to_reflected_method(cls::jclass, methodID::jmethodID, isStatic::jboolean, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ToReflectedMethod, jobject, (Ptr{JNIEnv}, jclass, jmethodID, jboolean,), penv, cls, methodID, isStatic)

get_superclass(sub::jclass, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetSuperclass, jclass, (Ptr{JNIEnv}, jclass,), penv, sub)

is_assignable_from(sub::jclass, sup::jclass, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].IsAssignableFrom, jboolean, (Ptr{JNIEnv}, jclass, jclass,), penv, sub, sup)

to_reflected_field(cls::jclass, fieldID::jfieldID, isStatic::jboolean, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ToReflectedField, jobject, (Ptr{JNIEnv}, jclass, jfieldID, jboolean,), penv, cls, fieldID, isStatic)

throw(obj::jthrowable, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].Throw, jint, (Ptr{JNIEnv}, jthrowable,), penv, obj)

throw_new(clazz::jclass, msg::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ThrowNew, jint, (Ptr{JNIEnv}, jclass, Cstring,), penv, clazz, msg)

exception_occurred(penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ExceptionOccurred, jthrowable, (Ptr{JNIEnv},), penv)

exception_describe(penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ExceptionDescribe, Nothing, (Ptr{JNIEnv},), penv)

exception_clear(penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ExceptionClear, Nothing, (Ptr{JNIEnv},), penv)

fatal_error(msg::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].FatalError, Nothing, (Ptr{JNIEnv}, Cstring,), penv, msg)

push_local_frame(capacity::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].PushLocalFrame, jint, (Ptr{JNIEnv}, jint,), penv, capacity)

pop_local_frame(result::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].PopLocalFrame, jobject, (Ptr{JNIEnv}, jobject,), penv, result)

new_global_ref(lobj::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewGlobalRef, jobject, (Ptr{JNIEnv}, jobject,), penv, lobj)

delete_global_ref(gref::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].DeleteGlobalRef, Nothing, (Ptr{JNIEnv}, jobject,), penv, gref)

delete_local_ref(obj::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].DeleteLocalRef, Nothing, (Ptr{JNIEnv}, jobject,), penv, obj)

is_same_object(obj1::jobject, obj2::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].IsSameObject, jboolean, (Ptr{JNIEnv}, jobject, jobject,), penv, obj1, obj2)

new_local_ref(ref::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewLocalRef, jobject, (Ptr{JNIEnv}, jobject,), penv, ref)

ensure_local_capacity(capacity::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].EnsureLocalCapacity, jint, (Ptr{JNIEnv}, jint,), penv, capacity)

alloc_object(clazz::jclass, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].AllocObject, jobject, (Ptr{JNIEnv}, jclass,), penv, clazz)

new_object_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewObjectA, jobject, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

get_object_class(obj::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetObjectClass, jclass, (Ptr{JNIEnv}, jobject,), penv, obj)

is_instance_of(obj::jobject, clazz::jclass, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].IsInstanceOf, jboolean, (Ptr{JNIEnv}, jobject, jclass,), penv, obj, clazz)

get_method_id(clazz::jclass, name::AbstractString, sig::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetMethodID, jmethodID, (Ptr{JNIEnv}, jclass, Cstring, Cstring,), penv, clazz, name, sig)

call_object_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallObjectMethodA, jobject, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_boolean_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallBooleanMethodA, jboolean, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_byte_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallByteMethodA, jbyte, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_char_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallCharMethodA, jchar, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_short_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallShortMethodA, jshort, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_int_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallIntMethodA, jint, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_long_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallLongMethodA, jlong, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_float_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallFloatMethodA, jfloat, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_double_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallDoubleMethodA, jdouble, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_void_method_a(obj::jobject, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallVoidMethodA, Nothing, (Ptr{JNIEnv}, jobject, jmethodID, Ptr{jvalue},), penv, obj, methodID, args)

call_nonvirtual_object_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualObjectMethodA, jobject, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

call_nonvirtual_boolean_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualBooleanMethodA, jboolean, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

call_nonvirtual_byte_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualByteMethodA, jbyte, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

call_nonvirtual_char_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualCharMethodA, jchar, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

call_nonvirtual_short_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualShortMethodA, jshort, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

call_nonvirtual_int_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualIntMethodA, jint, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

call_nonvirtual_long_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualLongMethodA, jlong, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

call_nonvirtual_float_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualFloatMethodA, jfloat, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

call_nonvirtual_double_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualDoubleMethodA, jdouble, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

call_nonvirtual_void_method_a(obj::jobject, clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallNonvirtualVoidMethodA, Nothing, (Ptr{JNIEnv}, jobject, jclass, jmethodID, Ptr{jvalue},), penv, obj, clazz, methodID, args)

get_field_id(clazz::jclass, name::AbstractString, sig::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetFieldID, jfieldID, (Ptr{JNIEnv}, jclass, Cstring, Cstring,), penv, clazz, name, sig)

get_object_field(obj::jobject, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetObjectField, jobject, (Ptr{JNIEnv}, jobject, jfieldID,), penv, obj, fieldID)

get_boolean_field(obj::jobject, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetBooleanField, jboolean, (Ptr{JNIEnv}, jobject, jfieldID,), penv, obj, fieldID)

get_byte_field(obj::jobject, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetByteField, jbyte, (Ptr{JNIEnv}, jobject, jfieldID,), penv, obj, fieldID)

get_char_field(obj::jobject, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetCharField, jchar, (Ptr{JNIEnv}, jobject, jfieldID,), penv, obj, fieldID)

get_short_field(obj::jobject, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetShortField, jshort, (Ptr{JNIEnv}, jobject, jfieldID,), penv, obj, fieldID)

get_int_field(obj::jobject, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetIntField, jint, (Ptr{JNIEnv}, jobject, jfieldID,), penv, obj, fieldID)

get_long_field(obj::jobject, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetLongField, jlong, (Ptr{JNIEnv}, jobject, jfieldID,), penv, obj, fieldID)

get_float_field(obj::jobject, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetFloatField, jfloat, (Ptr{JNIEnv}, jobject, jfieldID,), penv, obj, fieldID)

get_double_field(obj::jobject, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetDoubleField, jdouble, (Ptr{JNIEnv}, jobject, jfieldID,), penv, obj, fieldID)

set_object_field(obj::jobject, fieldID::jfieldID, val::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetObjectField, Nothing, (Ptr{JNIEnv}, jobject, jfieldID, jobject,), penv, obj, fieldID, val)

set_boolean_field(obj::jobject, fieldID::jfieldID, val::jboolean, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetBooleanField, Nothing, (Ptr{JNIEnv}, jobject, jfieldID, jboolean,), penv, obj, fieldID, val)

set_byte_field(obj::jobject, fieldID::jfieldID, val::jbyte, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetByteField, Nothing, (Ptr{JNIEnv}, jobject, jfieldID, jbyte,), penv, obj, fieldID, val)

set_char_field(obj::jobject, fieldID::jfieldID, val::jchar, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetCharField, Nothing, (Ptr{JNIEnv}, jobject, jfieldID, jchar,), penv, obj, fieldID, val)

set_short_field(obj::jobject, fieldID::jfieldID, val::jshort, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetShortField, Nothing, (Ptr{JNIEnv}, jobject, jfieldID, jshort,), penv, obj, fieldID, val)

set_int_field(obj::jobject, fieldID::jfieldID, val::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetIntField, Nothing, (Ptr{JNIEnv}, jobject, jfieldID, jint,), penv, obj, fieldID, val)

set_long_field(obj::jobject, fieldID::jfieldID, val::jlong, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetLongField, Nothing, (Ptr{JNIEnv}, jobject, jfieldID, jlong,), penv, obj, fieldID, val)

set_float_field(obj::jobject, fieldID::jfieldID, val::jfloat, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetFloatField, Nothing, (Ptr{JNIEnv}, jobject, jfieldID, jfloat,), penv, obj, fieldID, val)

set_double_field(obj::jobject, fieldID::jfieldID, val::jdouble, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetDoubleField, Nothing, (Ptr{JNIEnv}, jobject, jfieldID, jdouble,), penv, obj, fieldID, val)

get_static_method_id(clazz::jclass, name::AbstractString, sig::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticMethodID, jmethodID, (Ptr{JNIEnv}, jclass, Cstring, Cstring,), penv, clazz, name, sig)

call_static_object_method_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticObjectMethodA, jobject, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

call_static_boolean_method_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticBooleanMethodA, jboolean, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

call_static_byte_method_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticByteMethodA, jbyte, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

call_static_char_method_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticCharMethodA, jchar, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

call_static_short_method_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticShortMethodA, jshort, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

call_static_int_method_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticIntMethodA, jint, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

call_static_long_method_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticLongMethodA, jlong, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

call_static_float_method_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticFloatMethodA, jfloat, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

call_static_double_method_a(clazz::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticDoubleMethodA, jdouble, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, clazz, methodID, args)

call_static_void_method_a(cls::jclass, methodID::jmethodID, args::Array{jvalue,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].CallStaticVoidMethodA, Nothing, (Ptr{JNIEnv}, jclass, jmethodID, Ptr{jvalue},), penv, cls, methodID, args)

get_static_field_id(clazz::jclass, name::AbstractString, sig::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticFieldID, jfieldID, (Ptr{JNIEnv}, jclass, Cstring, Cstring,), penv, clazz, name, sig)

get_static_object_field(clazz::jclass, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticObjectField, jobject, (Ptr{JNIEnv}, jclass, jfieldID,), penv, clazz, fieldID)

get_static_boolean_field(clazz::jclass, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticBooleanField, jboolean, (Ptr{JNIEnv}, jclass, jfieldID,), penv, clazz, fieldID)

get_static_byte_field(clazz::jclass, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticByteField, jbyte, (Ptr{JNIEnv}, jclass, jfieldID,), penv, clazz, fieldID)

get_static_char_field(clazz::jclass, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticCharField, jchar, (Ptr{JNIEnv}, jclass, jfieldID,), penv, clazz, fieldID)

get_static_short_field(clazz::jclass, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticShortField, jshort, (Ptr{JNIEnv}, jclass, jfieldID,), penv, clazz, fieldID)

get_static_int_field(clazz::jclass, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticIntField, jint, (Ptr{JNIEnv}, jclass, jfieldID,), penv, clazz, fieldID)

get_static_long_field(clazz::jclass, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticLongField, jlong, (Ptr{JNIEnv}, jclass, jfieldID,), penv, clazz, fieldID)

get_static_float_field(clazz::jclass, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticFloatField, jfloat, (Ptr{JNIEnv}, jclass, jfieldID,), penv, clazz, fieldID)

get_static_double_field(clazz::jclass, fieldID::jfieldID, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStaticDoubleField, jdouble, (Ptr{JNIEnv}, jclass, jfieldID,), penv, clazz, fieldID)

set_static_object_field(clazz::jclass, fieldID::jfieldID, value::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetStaticObjectField, Nothing, (Ptr{JNIEnv}, jclass, jfieldID, jobject,), penv, clazz, fieldID, value)

set_static_boolean_field(clazz::jclass, fieldID::jfieldID, value::jboolean, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetStaticBooleanField, Nothing, (Ptr{JNIEnv}, jclass, jfieldID, jboolean,), penv, clazz, fieldID, value)

set_static_byte_field(clazz::jclass, fieldID::jfieldID, value::jbyte, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetStaticByteField, Nothing, (Ptr{JNIEnv}, jclass, jfieldID, jbyte,), penv, clazz, fieldID, value)

set_static_char_field(clazz::jclass, fieldID::jfieldID, value::jchar, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetStaticCharField, Nothing, (Ptr{JNIEnv}, jclass, jfieldID, jchar,), penv, clazz, fieldID, value)

set_static_short_field(clazz::jclass, fieldID::jfieldID, value::jshort, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetStaticShortField, Nothing, (Ptr{JNIEnv}, jclass, jfieldID, jshort,), penv, clazz, fieldID, value)

set_static_int_field(clazz::jclass, fieldID::jfieldID, value::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetStaticIntField, Nothing, (Ptr{JNIEnv}, jclass, jfieldID, jint,), penv, clazz, fieldID, value)

set_static_long_field(clazz::jclass, fieldID::jfieldID, value::jlong, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetStaticLongField, Nothing, (Ptr{JNIEnv}, jclass, jfieldID, jlong,), penv, clazz, fieldID, value)

set_static_float_field(clazz::jclass, fieldID::jfieldID, value::jfloat, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetStaticFloatField, Nothing, (Ptr{JNIEnv}, jclass, jfieldID, jfloat,), penv, clazz, fieldID, value)

set_static_double_field(clazz::jclass, fieldID::jfieldID, value::jdouble, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetStaticDoubleField, Nothing, (Ptr{JNIEnv}, jclass, jfieldID, jdouble,), penv, clazz, fieldID, value)

new_string(unicode::Array{jchar,1}, len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewString, jstring, (Ptr{JNIEnv}, Ptr{jchar}, jsize,), penv, unicode, len)

get_string_length(str::jstring, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStringLength, jsize, (Ptr{JNIEnv}, jstring,), penv, str)

get_string_chars(str::jstring, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStringChars, Ptr{jchar}, (Ptr{JNIEnv}, jstring, Ptr{jboolean},), penv, str, isCopy)

release_string_chars(str::jstring, chars::Array{jchar,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseStringChars, Nothing, (Ptr{JNIEnv}, jstring, Ptr{jchar},), penv, str, chars)

new_string_utf(utf::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewStringUTF, jstring, (Ptr{JNIEnv}, Cstring,), penv, utf)

get_string_utflength(str::jstring, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStringUTFLength, jsize, (Ptr{JNIEnv}, jstring,), penv, str)

get_string_utfchars(str::jstring, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStringUTFChars, Cstring, (Ptr{JNIEnv}, jstring, Ptr{jboolean},), penv, str, isCopy)

release_string_utfchars(str::jstring, chars::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseStringUTFChars, Nothing, (Ptr{JNIEnv}, jstring, Cstring,), penv, str, chars)

get_array_length(array::jarray, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetArrayLength, jsize, (Ptr{JNIEnv}, jarray,), penv, array)

new_object_array(len::Integer, clazz::jclass, init::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewObjectArray, jobjectArray, (Ptr{JNIEnv}, jsize, jclass, jobject,), penv, len, clazz, init)

get_object_array_element(array::jobjectArray, index::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetObjectArrayElement, jobject, (Ptr{JNIEnv}, jobjectArray, jsize,), penv, array, index)

set_object_array_element(array::jobjectArray, index::Integer, val::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetObjectArrayElement, Nothing, (Ptr{JNIEnv}, jobjectArray, jsize, jobject,), penv, array, index, val)

new_boolean_array(len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewBooleanArray, jbooleanArray, (Ptr{JNIEnv}, jsize,), penv, len)

new_byte_array(len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewByteArray, jbyteArray, (Ptr{JNIEnv}, jsize,), penv, len)

new_char_array(len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewCharArray, jcharArray, (Ptr{JNIEnv}, jsize,), penv, len)

new_short_array(len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewShortArray, jshortArray, (Ptr{JNIEnv}, jsize,), penv, len)

new_int_array(len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewIntArray, jintArray, (Ptr{JNIEnv}, jsize,), penv, len)

new_long_array(len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewLongArray, jlongArray, (Ptr{JNIEnv}, jsize,), penv, len)

new_float_array(len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewFloatArray, jfloatArray, (Ptr{JNIEnv}, jsize,), penv, len)

new_double_array(len::Integer, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewDoubleArray, jdoubleArray, (Ptr{JNIEnv}, jsize,), penv, len)

get_boolean_array_elements(array::jbooleanArray, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetBooleanArrayElements, Ptr{jboolean}, (Ptr{JNIEnv}, jbooleanArray, Ptr{jboolean},), penv, array, isCopy)

get_byte_array_elements(array::jbyteArray, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetByteArrayElements, Ptr{jbyte}, (Ptr{JNIEnv}, jbyteArray, Ptr{jboolean},), penv, array, isCopy)

get_char_array_elements(array::jcharArray, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetCharArrayElements, Ptr{jchar}, (Ptr{JNIEnv}, jcharArray, Ptr{jboolean},), penv, array, isCopy)

get_short_array_elements(array::jshortArray, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetShortArrayElements, Ptr{jshort}, (Ptr{JNIEnv}, jshortArray, Ptr{jboolean},), penv, array, isCopy)

get_int_array_elements(array::jintArray, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetIntArrayElements, Ptr{jint}, (Ptr{JNIEnv}, jintArray, Ptr{jboolean},), penv, array, isCopy)

get_long_array_elements(array::jlongArray, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetLongArrayElements, Ptr{jlong}, (Ptr{JNIEnv}, jlongArray, Ptr{jboolean},), penv, array, isCopy)

get_float_array_elements(array::jfloatArray, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetFloatArrayElements, Ptr{jfloat}, (Ptr{JNIEnv}, jfloatArray, Ptr{jboolean},), penv, array, isCopy)

get_double_array_elements(array::jdoubleArray, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetDoubleArrayElements, Ptr{jdouble}, (Ptr{JNIEnv}, jdoubleArray, Ptr{jboolean},), penv, array, isCopy)

release_boolean_array_elements(array::jbooleanArray, elems::Ptr{jboolean}, mode::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseBooleanArrayElements, Nothing, (Ptr{JNIEnv}, jbooleanArray, Ptr{jboolean}, jint,), penv, array, elems, mode)

release_byte_array_elements(array::jbyteArray, elems::Ptr{jbyte}, mode::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseByteArrayElements, Nothing, (Ptr{JNIEnv}, jbyteArray, Ptr{jbyte}, jint,), penv, array, elems, mode)

release_char_array_elements(array::jcharArray, elems::Ptr{jchar}, mode::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseCharArrayElements, Nothing, (Ptr{JNIEnv}, jcharArray, Ptr{jchar}, jint,), penv, array, elems, mode)

release_short_array_elements(array::jshortArray, elems::Ptr{jshort}, mode::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseShortArrayElements, Nothing, (Ptr{JNIEnv}, jshortArray, Ptr{jshort}, jint,), penv, array, elems, mode)

release_int_array_elements(array::jintArray, elems::Ptr{jint}, mode::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseIntArrayElements, Nothing, (Ptr{JNIEnv}, jintArray, Ptr{jint}, jint,), penv, array, elems, mode)

release_long_array_elements(array::jlongArray, elems::Ptr{jlong}, mode::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseLongArrayElements, Nothing, (Ptr{JNIEnv}, jlongArray, Ptr{jlong}, jint,), penv, array, elems, mode)

release_float_array_elements(array::jfloatArray, elems::Ptr{jfloat}, mode::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseFloatArrayElements, Nothing, (Ptr{JNIEnv}, jfloatArray, Ptr{jfloat}, jint,), penv, array, elems, mode)

release_double_array_elements(array::jdoubleArray, elems::Ptr{jdouble}, mode::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseDoubleArrayElements, Nothing, (Ptr{JNIEnv}, jdoubleArray, Ptr{jdouble}, jint,), penv, array, elems, mode)

get_boolean_array_region(array::jbooleanArray, start::Integer, l::Integer, buf::Array{jboolean,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetBooleanArrayRegion, Nothing, (Ptr{JNIEnv}, jbooleanArray, jsize, jsize, Ptr{jboolean},), penv, array, start, l, buf)

get_byte_array_region(array::jbyteArray, start::Integer, len::Integer, buf::Array{jbyte,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetByteArrayRegion, Nothing, (Ptr{JNIEnv}, jbyteArray, jsize, jsize, Ptr{jbyte},), penv, array, start, len, buf)

get_char_array_region(array::jcharArray, start::Integer, len::Integer, buf::Array{jchar,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetCharArrayRegion, Nothing, (Ptr{JNIEnv}, jcharArray, jsize, jsize, Ptr{jchar},), penv, array, start, len, buf)

get_short_array_region(array::jshortArray, start::Integer, len::Integer, buf::Array{jshort,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetShortArrayRegion, Nothing, (Ptr{JNIEnv}, jshortArray, jsize, jsize, Ptr{jshort},), penv, array, start, len, buf)

get_int_array_region(array::jintArray, start::Integer, len::Integer, buf::Array{jint,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetIntArrayRegion, Nothing, (Ptr{JNIEnv}, jintArray, jsize, jsize, Ptr{jint},), penv, array, start, len, buf)

get_long_array_region(array::jlongArray, start::Integer, len::Integer, buf::Array{jlong,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetLongArrayRegion, Nothing, (Ptr{JNIEnv}, jlongArray, jsize, jsize, Ptr{jlong},), penv, array, start, len, buf)

get_float_array_region(array::jfloatArray, start::Integer, len::Integer, buf::Array{jfloat,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetFloatArrayRegion, Nothing, (Ptr{JNIEnv}, jfloatArray, jsize, jsize, Ptr{jfloat},), penv, array, start, len, buf)

get_double_array_region(array::jdoubleArray, start::Integer, len::Integer, buf::Array{jdouble,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetDoubleArrayRegion, Nothing, (Ptr{JNIEnv}, jdoubleArray, jsize, jsize, Ptr{jdouble},), penv, array, start, len, buf)

set_boolean_array_region(array::jbooleanArray, start::Integer, l::Integer, buf::Array{jboolean,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetBooleanArrayRegion, Nothing, (Ptr{JNIEnv}, jbooleanArray, jsize, jsize, Ptr{jboolean},), penv, array, start, l, buf)

set_byte_array_region(array::jbyteArray, start::Integer, len::Integer, buf::Array{jbyte,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetByteArrayRegion, Nothing, (Ptr{JNIEnv}, jbyteArray, jsize, jsize, Ptr{jbyte},), penv, array, start, len, buf)

set_char_array_region(array::jcharArray, start::Integer, len::Integer, buf::Array{jchar,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetCharArrayRegion, Nothing, (Ptr{JNIEnv}, jcharArray, jsize, jsize, Ptr{jchar},), penv, array, start, len, buf)

set_short_array_region(array::jshortArray, start::Integer, len::Integer, buf::Array{jshort,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetShortArrayRegion, Nothing, (Ptr{JNIEnv}, jshortArray, jsize, jsize, Ptr{jshort},), penv, array, start, len, buf)

set_int_array_region(array::jintArray, start::Integer, len::Integer, buf::Array{jint,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetIntArrayRegion, Nothing, (Ptr{JNIEnv}, jintArray, jsize, jsize, Ptr{jint},), penv, array, start, len, buf)

set_long_array_region(array::jlongArray, start::Integer, len::Integer, buf::Array{jlong,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetLongArrayRegion, Nothing, (Ptr{JNIEnv}, jlongArray, jsize, jsize, Ptr{jlong},), penv, array, start, len, buf)

set_float_array_region(array::jfloatArray, start::Integer, len::Integer, buf::Array{jfloat,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetFloatArrayRegion, Nothing, (Ptr{JNIEnv}, jfloatArray, jsize, jsize, Ptr{jfloat},), penv, array, start, len, buf)

set_double_array_region(array::jdoubleArray, start::Integer, len::Integer, buf::Array{jdouble,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].SetDoubleArrayRegion, Nothing, (Ptr{JNIEnv}, jdoubleArray, jsize, jsize, Ptr{jdouble},), penv, array, start, len, buf)

register_natives(clazz::jclass, methods::Array{JNINativeMethod,1}, nMethods::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].RegisterNatives, jint, (Ptr{JNIEnv}, jclass, Ptr{JNINativeMethod}, jint,), penv, clazz, methods, nMethods)

unregister_natives(clazz::jclass, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].UnregisterNatives, jint, (Ptr{JNIEnv}, jclass,), penv, clazz)

monitor_enter(obj::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].MonitorEnter, jint, (Ptr{JNIEnv}, jobject,), penv, obj)

monitor_exit(obj::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].MonitorExit, jint, (Ptr{JNIEnv}, jobject,), penv, obj)

get_java_vm(vm::Array{JavaVM,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetJavaVM, jint, (Ptr{JNIEnv}, Array{JavaVM,1},), penv, vm)

get_string_region(str::jstring, start::Integer, len::Integer, buf::Array{jchar,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStringRegion, Nothing, (Ptr{JNIEnv}, jstring, jsize, jsize, Ptr{jchar},), penv, str, start, len, buf)

get_string_utfregion(str::jstring, start::Integer, len::Integer, buf::AbstractString, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStringUTFRegion, Nothing, (Ptr{JNIEnv}, jstring, jsize, jsize, Cstring,), penv, str, start, len, buf)

get_primitive_array_critical(array::jarray, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetPrimitiveArrayCritical, Ptr{Nothing}, (Ptr{JNIEnv}, jarray, Ptr{jboolean},), penv, array, isCopy)

release_primitive_array_critical(array::jarray, carray::Ptr{Nothing}, mode::jint, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleasePrimitiveArrayCritical, Nothing, (Ptr{JNIEnv}, jarray, Ptr{Nothing}, jint,), penv, array, carray, mode)

get_string_critical(string::jstring, isCopy::Ptr{jboolean}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetStringCritical, Ptr{jchar}, (Ptr{JNIEnv}, jstring, Ptr{jboolean},), penv, string, isCopy)

release_string_critical(string::jstring, cstring::Array{jchar,1}, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ReleaseStringCritical, Nothing, (Ptr{JNIEnv}, jstring, Ptr{jchar},), penv, string, cstring)

new_weak_global_ref(obj::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewWeakGlobalRef, jweak, (Ptr{JNIEnv}, jobject,), penv, obj)

delete_weak_global_ref(ref::jweak, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].DeleteWeakGlobalRef, Nothing, (Ptr{JNIEnv}, jweak,), penv, ref)

exception_check(penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].ExceptionCheck, jboolean, (Ptr{JNIEnv},), penv)

new_direct_byte_buffer(address::Ptr{Nothing}, capacity::jlong, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].NewDirectByteBuffer, jobject, (Ptr{JNIEnv}, Ptr{Nothing}, jlong,), penv, address, capacity)

get_direct_buffer_address(buf::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetDirectBufferAddress, Ptr{Nothing}, (Ptr{JNIEnv}, jobject,), penv, buf)

get_direct_buffer_capacity(buf::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetDirectBufferCapacity, jlong, (Ptr{JNIEnv}, jobject,), penv, buf)

get_object_ref_type(obj::jobject, penv=jnienvptrs[]) =
  ccall(jninativeinterfaceref[].GetObjectRefType, jobjectRefType, (Ptr{JNIEnv}, jobject,), penv, obj)

# === Above Generated by make_jni2.jl ===


text:
using Base.Test
using MachineLearningMetrics

y_true = [1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0]
y_pred = [1, 1, 0, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1]

# total_population
@test total_population(y_true, y_pred) == 16

# true_positive
@test true_positive(y_true, y_pred) == 4

# true_negative
@test true_negative(y_true, y_pred) == 4

# false_positive
@test false_positive(y_true, y_pred) == 4

# false_negative
@test false_negative(y_true, y_pred) == 4

# condition_positive
@test condition_positive(y_true, y_pred) == 8

# condition_negative
@test condition_negative(y_true, y_pred) == 8

# predicted_condition_positive
@test predicted_condition_positive(y_true, y_pred) == 8

# predicted_condition_negative
@test predicted_condition_negative(y_true, y_pred) == 8

# accuracy_score
# accuracy
@test accuracy_score(y_true, y_pred) == accuracy(y_true, y_pred) == 0.5

# prevalence
@test accuracy(y_true, y_pred) == 0.5

# positive_predictive_value
# precision
@test positive_predictive_value(y_true, y_pred) == precision(y_true, y_pred) == 0.5

# false_discovery_rate
@test false_discovery_rate(y_true, y_pred) == 0.5

# negative_predictive_value
@test negative_predictive_value(y_true, y_pred) == 0.5

# false_omission_rate
@test false_omission_rate(y_true, y_pred) == 0.5

# true_positive_rate
# sensitivity
# recall
@test positive_predictive_value(y_true, y_pred) == precision(y_true, y_pred) == 0.5

# false_positive_rate
@test false_positive_rate(y_true, y_pred) == 0.5

# false_negative_rate
@test false_negative_rate(y_true, y_pred) == 0.5

# true_negative_rate
# specificity
@test true_negative_rate(y_true, y_pred) == specificity(y_true, y_pred) == 0.5

# positive_likelihood_ratio
@test positive_likelihood_ratio(y_true, y_pred) == 1.0

# negative_likelihood_ratio
@test negative_likelihood_ratio(y_true, y_pred) == 1.0

# diagnostic_odds_ratio
@test diagnostic_odds_ratio(y_true, y_pred) == 1.0

# f1_score
@test f1_score(y_true, y_pred) == 0.5

# matthews_corrcoef
@test matthews_corrcoef(y_true, y_pred) == 0.0


text:
"""
Format.jl

This package provides various functions to provide formatted output,
either in a fashion similar to C printf or Python format strings.

## Python-style Types and Functions

#### Types to Represent Formats

This package has two types ``FormatSpec`` and ``FormatExpr`` to represent a format specification.

In particular, ``FormatSpec`` is used to capture the specification of a single entry. One can compile a format specification string into a ``FormatSpec`` instance as

```julia
fspec = FormatSpec("d")
fspec = FormatSpec("<8.4f")
```
Please refer to [Python's format specification language](http://docs.python.org/2/library/string.html#formatspec) for details.

``FormatExpr`` captures a formatting expression that may involve multiple items. One can compile a formatting string into a ``FormatExpr`` instance as

```julia
fe = FormatExpr("{1} + {2}")
fe = FormatExpr("{1:d} + {2:08.4e} + {3|>abs2}")
```
Please refer to [Python's format string syntax](http://docs.python.org/2/library/string.html#format-string-syntax) for details.

**Note:** If the same format is going to be applied for multiple times. It is more efficient to first compile it.

#### Formatted Printing

One can use ``printfmt`` and ``printfmtln`` for formatted printing:

- **printfmt**(io, fe, args...)

- **printfmt**(fe, args...)

    Print given arguments using given format ``fe``. Here ``fe`` can be a formatting string, an instance of ``FormatSpec`` or ``FormatExpr``.

    **Examples**

    ```julia
    printfmt("{1:>4s} + {2:.2f}", "abc", 12) # --> print(" abc + 12.00")
    printfmt("{} = {:#04x}", "abc", 12) # --> print("abc = 0x0c")

    fs = FormatSpec("#04x")
    printfmt(fs, 12)   # --> print("0x0c")

    fe = FormatExpr("{} = {:#04x}")
    printfmt(fe, "abc", 12)   # --> print("abc = 0x0c")
    ```

    **Notes**

    If the first argument is a string, it will be first compiled into a ``FormatExpr``, which implies that you can not use specification-only string in the first argument.

    ```julia
    printfmt("{1:d}", 10)   # OK, "{1:d}" can be compiled into a FormatExpr instance
    printfmt("d", 10)       # Error, "d" can not be compiled into a FormatExpr instance
                            # such a string to specify a format specification for single argument

    printfmt(FormatSpec("d"), 10)  # OK
    printfmt(FormatExpr("{1:d}", 10)) # OK
    ```

- **printfmtln**(io, fe, args...)

- **printfmtln**(fe, args...)

    Similar to ``printfmt`` except that this function print a newline at the end.

#### Formatted String

One can use ``pyfmt`` to format a single value into a string, or ``format`` to format one to multiple arguments into a string using an format expression.

- **pyfmt**(fspec, a)

    Format a single value using a format specification given by ``fspec``, where ``fspec`` can be either a string or an instance of ``FormatSpec``.

- **format**(fe, args...)

    Format arguments using a format expression given by ``fe``, where ``fe`` can be either a string or an instance of ``FormatSpec``.

#### Difference from Python's Format

At this point, this package implements a subset of Python's formatting language (with slight modification). Here is a summary of the differences:

- ``g`` and ``G`` for floating point formatting have not been supported yet. Please use ``f``, ``e``, or ``E`` instead.

- The package currently provides default alignment, left alignment ``<`` and right alignment ``>``. Other form of alignment such as centered alignment ``^`` has not been supported yet.

- In terms of argument specification, it supports natural ordering (e.g. ``{} + {}``), explicit position (e.g. ``{1} + {2}``). It hasn't supported named arguments or fields extraction yet. Note that mixing these two modes is not allowed (e.g. ``{1} + {}``).

- The package provides support for filtering (for explicitly positioned arguments), such as ``{1|>lowercase}`` by allowing one to embed the ``|>`` operator, which the Python counter part does not support.

## C-style functions

The c-style part of this package aims to get around the limitation that
`@sprintf` has to take a literal string argument.
The core part is basically a c-style print formatter using the standard
`@sprintf` macro.
It also adds functionalities such as commas separator (thousands), parenthesis for negatives,
stripping trailing zeros, and mixed fractions.

### Usage and Implementation

The idea here is that the package compiles a function only once for each unique
format string within the `Format.*` name space, so repeated use is faster.
Unrelated parts of a session using the same format string would reuse the same
function, avoiding redundant compilation. To avoid the proliferation of
functions, we limit the usage to only 1 argument. Practical consideration
would suggest that only dozens of functions would be created in a session, which
seems manageable.

Usage
```julia
using Format

fmt = "%10.3f"
s = cfmt( fmt, 3.14159 ) # usage 1. Quite performant. Easiest to switch to.

fmtrfunc = generate_formatter( fmt ) # usage 2. This bypass repeated lookup of cached function. Most performant.
s = fmtrfunc( 3.14159 )

s = format( 3.14159, precision=3 ) # usage 3. Most flexible, with some non-printf options. Least performant.
```
### Speed

`cfmt`: Speed penalty is about 20% for floating point and 30% for integers.

If the formatter is stored and used instead (see the example using `generate_formatter` above),
the speed penalty reduces to 10% for floating point and 15% for integers.

### Commas

This package also supplements the lack of thousand separator e.g. `"%'d"`, `"%'f"`, `"%'s"`.

Note: `"%'s"` behavior is that for small enough floating point (but not too small),
thousand separator would be used. If the number needs to be represented by `"%e"`, no
separator is used.

### Flexible `format` function

This package contains a run-time number formatter `format` function, which goes beyond
the standard `sprintf` functionality.

An example:
```julia
s = format( 1234, commas=true ) # 1,234
s = format( -1234, commas=true, parens=true ) # (1,234)
```

The keyword arguments are (Bold keywards are not printf standard)

* width. Integer. Try to fit the output into this many characters. May not be successful.
   Sacrifice space first, then commas.
* precision. Integer. How many decimal places.
* leftjustified. Boolean
* zeropadding. Boolean
* commas. Boolean. Thousands-group separator.
* signed. Boolean. Always show +/- sign?
* positivespace. Boolean. Prepend an extra space for positive numbers? (so they align nicely with negative numbers)
* **parens**. Boolean. Use parenthesis instead of "-". e.g. `(1.01)` instead of `-1.01`. Useful in finance. Note that
  you cannot use `signed` and `parens` option at the same time.
* **stripzeros**. Boolean. Strip trailing '0' to the right of the decimal (and to the left of 'e', if any ).
   * It may strip the decimal point itself if all trailing places are zeros.
   * This is true by default if precision is not given, and vice versa.
* alternative. Boolean. See `#` alternative form explanation in standard printf documentation
* conversion. length=1 string. Default is type dependent. It can be one of `aAeEfFoxX`. See standard
  printf documentation.
* **mixedfraction**. Boolean. If the number is rational, format it in mixed fraction e.g. `1_1/2` instead of `3/2`
* **mixedfractionsep**. Default `_`
* **fractionsep**. Default `/`
* **fractionwidth**. Integer. Try to pad zeros to the numerator until the fractional part has this width
* **tryden**. Integer. Try to use this denominator instead of a smaller one. No-op if it'd lose precision.
* **suffix**. String. This strings will be appended to the output. Useful for units/%
* **autoscale**. Symbol, default `:none`. It could be `:metric`, `:binary`, or `:finance`.
    * `:metric` implements common SI symbols for large and small numbers e.g. `M`, `k`, `μ`, `n`
    * `:binary` implements common ISQ symbols for large numbers e.g. `Ti`, `Gi`, `Mi`, `Ki`
    * `:finance` implements common finance/news symbols for large numbers e.g. `b` (billion), `m` (millions)
"""
module Format

import Base.show

_stdout() = stdout
_codeunits(s) = Vector{UInt8}(codeunits(s))
m_eval(expr) = Core.eval(@__MODULE__, expr)

export FormatSpec, FormatExpr, printfmt, printfmtln, format, generate_formatter
export pyfmt, cfmt, fmt
export fmt_default, fmt_default!, reset!, default_spec, default_spec!

# Deal with mess from #16058
# Later, use Strs package!
isdefined(Main, :ASCIIStr) || (const ASCIIStr = String)
isdefined(Main, :UTF8Str)  || (const UTF8Str = String)

include("cformat.jl" )
include("fmtspec.jl")
include("fmtcore.jl")
include("formatexpr.jl")
include("fmt.jl")

end # module Format


text:
using MLDataStructs
using Test

@testset "MLDataStructs.jl" begin
    # Write your own tests here.
end


text:
# --------------------------------------------------------------------------- #

# Setting an ip model of SPP
function setSPP(solverSelected, C, A)
  m, n = size(A)
  ip = Model(solver=solverSelected)
  @variable(ip, x[1:n], Bin)
  @objective(ip, Max, dot(C, x))
  @constraint(ip , cte[i=1:m], sum(A[i,j] * x[j] for j=1:n) <= 1)
  return ip, x
end


text:
@doc doc"""
Ising model with energy $E = -\sum_{ij} J_{ij} \sigma_i \sigma_j$,
where $\sigma_i$ takes value of 1 (up spin) or -1 (down spin).
"""
mutable struct Ising <: Model
    lat :: Lattice
    spins :: Matrix{Int}
    rng :: Random.MersenneTwister

    function Ising(lat::Lattice, rng::Random.AbstractRNG)
        model = new()
        model.lat = lat
        model.rng = rng
        model.spins = rand(model.rng, [1,-1], 1, numsites(lat))
        return model
    end
end

Ising(lat::Lattice) = Ising(lat, Random.seed!(Random.MersenneTwister(0)))
Ising(lat::Lattice, seed) = Ising(lat, Random.seed!(Random.MersenneTwister(0), seed...))

@doc doc"""
    Ising(param)

Generates `Ising` using `param["Lattice"]` and `param["Seed"]` (if defined).
Each spin will be initialized randomly and independently.
"""
function Ising(param::Parameter)
    lat = generatelattice(param)
    if "Seed" in keys(param)
        return Ising(lat, param["Seed"])
    else
        return Ising(lat)
    end
end

include("update.jl")
include("estimator.jl")
include("postproc.jl")


text:
using PyTensorFlow
using PyCall
using LinearAlgebra
using PyPlot
using Random
Random.seed!(233)

if Sys.islinux()
py"""
import tensorflow as tf
libLaplacian = tf.load_op_library('../Ops/Laplacian/build/libLaplacian.so')
@tf.custom_gradient
def laplacian_op(coef,func,h,rhograv):
    p = libLaplacian.laplacian(coef,func,h,rhograv)
    def grad(dy):
        return libLaplacian.laplacian_grad(dy, coef, func, h, rhograv)
    return p, grad
"""
elseif Sys.isapple()
py"""
import tensorflow as tf
libPoissonOp = tf.load_op_library('../Ops/Laplacian/build/libLaplacian.dylib')
@tf.custom_gradient
def laplacian_op(coef,func,h,rhograv):
    p = libLaplacian.laplacian(coef,func,h,rhograv)
    def grad(dy):
        return libLaplacian.laplacian_grad(dy, coef, func, h, rhograv)
    return p, grad
"""
elseif Sys.iswindows()
py"""
import tensorflow as tf
libPoissonOp = tf.load_op_library('../Ops/Laplacian/build/libLaplacian.dll')
@tf.custom_gradient
def laplacian_op(coef,func,h,rhograv):
    p = libLaplacian.laplacian(coef,func,h,rhograv)
    def grad(dy):
        return libLaplacian.laplacian_grad(dy, coef, func, h, rhograv)
    return p, grad
"""
end

laplacian_op = py"laplacian_op"

if Sys.islinux()
py"""
import tensorflow as tf
libUpwlapOp = tf.load_op_library('../Ops/Upwlap/build/libUpwlapOp.so')
@tf.custom_gradient
def upwlap_op(perm,mobi,func,h,rhograv):
    out = libUpwlapOp.upwlap_op(perm,mobi,func,h,rhograv)
    def grad(dy):
        return libUpwlapOp.upwlap_op_grad(dy, out, perm,mobi,func,h,rhograv)
    return out, grad
"""
elseif Sys.isapple()
py"""
import tensorflow as tf
libUpwlapOp = tf.load_op_library('../Ops/Upwlap/build/libUpwlapOp.dylib')
@tf.custom_gradient
def upwlap_op(perm,mobi,func,h,rhograv):
    out = libUpwlapOp.upwlap_op(perm,mobi,func,h,rhograv)
    def grad(dy):
        return libUpwlapOp.upwlap_op_grad(dy, out, perm,mobi,func,h,rhograv)
    return out, grad
"""
elseif Sys.iswindows()
py"""
import tensorflow as tf
libUpwlapOp = tf.load_op_library('./Ops/Upwlap/build/libUpwlapOp.dll')
@tf.custom_gradient
def upwlap_op(perm,mobi,func,h,rhograv):
    out = libUpwlapOp.upwlap_op(perm,mobi,func,h,rhograv)
    def grad(dy):
        return libUpwlapOp.upwlap_op_grad(dy, out, perm,mobi,func,h,rhograv)
    return out, grad
"""
end

upwlap_op = py"upwlap_op"


text:
include("imports.jl")

@testset ExtendedTestSet "OneOf" begin
    item = ArrayItem(rand(10, 10))
    tfm = OneOf([MapElem(x -> x + 1), MapElem(x -> x - 1)])
    @test apply(tfm, item; randstate = (1, nothing)).data ≈ item.data .+ 1
    @test apply(tfm, item; randstate = (2, nothing)).data ≈ item.data .- 1
    testapply(tfm, ArrayItem)
    testapply!(tfm, ArrayItem)
end

@testset ExtendedTestSet "Maybe" begin
    item = ArrayItem(rand(10, 10))
    tfm = Maybe(MapElem(x -> x + 1))
    @test apply(tfm, item; randstate = (1, nothing)).data ≈ item.data .+ 1
    @test apply(tfm, item; randstate = (2, nothing)).data ≈ item.data
    testapply(tfm, ArrayItem)
end

@testset ExtendedTestSet "OneOfProjective" begin
    item = Image(rand(10, 10))
    tfm = Rotate(10) |> CenterCrop((5, 5))
    oneof = Maybe(tfm)
    @test oneof isa DataAugmentation.OneOfProjective
    @test_nowarn DataAugmentation.getprojection(oneof, getbounds(item); randstate = getrandstate(oneof))
    @test apply(oneof, item; randstate = (1, getrandstate(tfm))) isa Image
    @test apply(oneof, item; randstate = (2, nothing)) isa Image
end

@testset ExtendedTestSet "Sequence(oneofs...)" begin
    item = ArrayItem(rand(10, 10))
    tfm = Maybe(MapElem(x -> x + 1)) |> Maybe(MapElem(x -> x + 1))
    testapply(tfm, ArrayItem)
    testapply!(tfm, ArrayItem)
end

@testset ExtendedTestSet "Sequence" begin
    tfm = Maybe(AdjustBrightness(0.1), .75) |> Maybe(AdjustContrast(0.1), .75)
    item = DataAugmentation.Image(rand(RGB{Float32}, 10, 10))
    testapply(tfm, Image)
    testapply!(tfm, Image)
end