text:
using LangrethsRules
using Test

function integrate(x::AbstractVector, y::AbstractVector)
    if isone(length(y))
        return zero(first(y))
    end
    @inbounds retval = (x[2] - x[1]) * (y[1] + y[2])
    @inbounds @fastmath @simd for i in 2:(length(y)-1)
        retval += (x[i+1] - x[i]) * (y[i] + y[i+1])
    end
    return 1 // 2 * retval
end

function compute(A, B, ts)
    l = zero(lesser(A))
    for i in 1:size(l, 1)
        for j in 1:i
            l[i,j] += integrate(ts, greater(A)[i, 1:i] .* lesser(B)[1:i, j])
            l[i,j] -= integrate(ts, lesser(A)[i, 1:j] .* greater(B)[1:j, j])
            l[i,j] -= integrate(ts[j:i], lesser(A)[i, j:i] .* lesser(B)[j:i, j])

            l[j,i] = -conj(l[i,j])
        end
        l[i,i] = eltype(l) <: Real ? 0.0 : im * imag(l[i,i])
    end

    g = zero(greater(A))
    for i in 1:size(g, 1)
        for j in 1:i
            g[i,j] -= integrate(ts, lesser(A)[i, 1:i] .* greater(B)[1:i, j])
            g[i,j] += integrate(ts, greater(A)[i, 1:j] .* lesser(B)[1:j, j])
            g[i,j] += integrate(ts[j:i], greater(A)[i, j:i] .* greater(B)[j:i, j])

            g[j,i] = -conj(g[i,j])
        end
        g[i,i] = eltype(g) <: Real ? 0.0 : im * imag(g[i,i])
    end

    TimeOrderedGF(l, g, ts)
end

N = 100

# Suppose a non-equidistant time-grid
ts = sort(N*rand(N));

# And 2 time-ordered GFs defined on that grid
a = let
    x = rand(ComplexF64,N,N)
    y = rand(ComplexF64,N,N)
    TimeOrderedGF(x - x', y - y', ts) # Skew-symmetric L and G components
end
b = let
    x = rand(ComplexF64,N,N)
    y = rand(ComplexF64,N,N)
    TimeOrderedGF(x - x', y - y', ts) # Skew-symmetric L and G components
end

@test greater(a ⋆ b) ≈ greater(compute(a, b, ts))
@test lesser(a ⋆ b) ≈ lesser(compute(a, b, ts))
@test retarded(a ⋆ b) ≈ adjoint(advanced(a ⋆ b))
# @test retarded(a ⋆ b ⋆ a) ≈ retarded(a ⋆ (b ⋆ a))


text:
import LsqFit

# #############################################################################
#Moving average filtering

#Apply filter
function movavgfilt{T}(y::Array{T,1},N::Int64)
    #Usage:
    #   ys = movavgfilt(y,N)
    #       y       Signal to be smoothed
    #       N       Size of window
    #       ys      Smoothed output

    #Initiate
    ys::Array{Float64,1} = Array(Float64,length(y))

    #Compute mean
    for i=1:length(y)
        s::Int64 = minimum([floor((N-1)/2),i-1,length(y)-i])
        ys[i] = sum(y[i-s:i+s])/(2*s+1)
    end

    return ys
end

# #############################################################################
#Savitzsky-Golay filtering

#Get coefficients
function sgolaycoeff{T<:Real}(z::Array{T,1},deg::Int64,diff::Int64=0)
    #Usage:
    #   C = sgolaycoeff(z,deg)
    #   C = sgolaycoeff(z,deg,diff)
    #       z       Fitting window (eg: [-5:5])
    #       deg     Degree of polynomial for fit
    #       diff    Number of differentiations (must be <= deg, default 0 for
    #               smoothing)
    #       C       Matrix for transforming window values to polynomial coefficients,
    #               where (i,j)th entry is the multiplier to the jth value for the
    #               (i-1) order term in the fitting polynomial.

    #Check input
    assert(diff <= deg)

    #Solve normal equations
    J::Array{T,2} = broadcast(^,z,transpose(0:deg))
    Jt::Array{T,2} = transpose(J)
    C::Array{Float64,2} = inv(Jt*J)*Jt

    #Get coefficients for derivatives
    for i = 1:diff
        C = broadcast(*,[1:size(C,1)-1],C[2:end,:])
    end

    return C
end

#Apply filter
function sgolayfilt{T<:Real}(y::Array{T,1},N::Int64;deg::Int64=2,diff::Int64=0,multi::Int64=1)
    #Usage:
    #   ys = sgolayfilt(y,N)
    #   ys = sgolayfilt(x,y,N)
    #   ys = sgolayfilt(..., deg=deg)
    #   ys = sgolayfilt(..., diff=diff)
    #   ys = sgolayfilt(..., multi=multi)
    #       y       Signal to be smoothed
    #       N       Size of window
    #       x       Distribution of signal events (default [1:length(y)] for uniform)
    #       deg     Degree of polynomial for fit
    #       diff    Number of differentiations (must be <= deg, default 0 for
    #               smoothing)
    #       multi   Number of times to filter
    #       ys      Smoothed signal

    #Check input
    if diff > deg
        error("Differential order cannot exceed polynomial fit degree.")
    end
    m::Int64 = div(N-1,2)
    if 2*m+1 <= deg
        error("Fitting window N = $N is too small for polynomial fit of degree $deg.")
    end
    if multi < 1
        error("Number of times to filter must be at least 1.")
    end

    #Initiate
    C::Array{Float64,2} = sgolaycoeff([-m:m],deg,diff)
    ys::Array{Float64,1} = Array(Float64,length(y))

    #Starting
    a::Array{Float64,1} = C*y[1:2m+1]
    for i=1:m
        ys[i] = sum(a.*(i-m-1).^[0:length(a)-1])
    end

    #Middle
    C1::Array{Float64,1} = C[1,:][:]
    for i=m+1:length(y)-m
        ys[i] = dot(C1,y[i-m:i+m])
    end

    #Ending
    a = C*y[end-2m:end]
    for i=length(y)-m+1:length(y)
        ys[i] = sum(a.*(i-length(y)+m).^[0:length(a)-1])
    end

    if multi == 1
        return ys
    else
        return sgolayfilt(ys,N,deg=deg,diff=0,multi=multi-1)
    end
end

function sgolayfilt{T<:Real}(x::Array{T,1},y::Array{T,1},N::Int64;deg::Int64=2,diff::Int64=0)
    #Check input
    if diff > deg
        error("Differential order cannot exceed polynomial fit degree.")
    end
    m::Int64 = div(N-1,2)
    if 2*m+1 <= deg
        error("Fitting window N = $N is too small for polynomial fit of degree $deg.")
    end
    if multi < 1
        error("Number of times to filter must be at least 1.")
    end

    #Initiate
    ys::Array{Float64,1} = Array(Float64,length(y))

    #Starting
    C::Array{Float64,2} = sgolaycoeff(x[1:2m+1],deg,diff)
    a::Array{Float64,1} = C*y[1:2m+1]
    for i=1:m
        ys[i] = sum(a.*x[i].^[0:length(a)-1])
    end

    #Middle
    for i=m+1:length(y)-m
        C = sgolaycoeff(x[i-m:i+m],deg,diff)
        a = C*y[i-m:i+m]
        ys[i] = sum(a.*x[i].^[0:length(a)-1])
    end

    #Ending
    C = sgolaycoeff(x[end-2*m:end],deg,diff)
    a = C*y[end-2*m:end]
    for i=length(y)-m+1:length(y)
        ys[i] = sum(a.*x[i].^[0:length(a)-1])
    end

    if multi == 1
        return ys
    else
        return sgolayfilt(ys,N,deg=deg,diff=0,multi=multi-1)
    end
end

# #############################################################################
#Robust Local Weighted Scatter-plot Smooth
function robust_lowess{T}(y::Array{T,1},N::Int64;deg::Int64=1,robust::Bool=true)
    #Initiate
    ys::Array{Float64,1} = Array(Float64,length(y))
    w::Array{Float64,1} = Array(Float64,N)
    s::Int64 = floor((N-1)/2)

    #Modelling function
    f(x::Array,p::Array{Float64,1}) = sum(broadcast(^,x,transpose([0:length(p)-1]))*p,2)
    f(x::Int64,p::Array{Float64,1}) = sum(p.*x.^[0:length(p)-1])

    #Conduct regression
    for i=1:length(y)
        #Get regression region
        if s > i-1
            range::Array{Int64,1} = [1:N]
        elseif s > length(y)-i
            range = [(length(y)-N+1):length(y)]
        else
            range = [(i-s):(i+s)]
        end

        #Compute local weights
        w = (1-(abs(range-i)/maximum(abs(range-i))).^3).^3

        #Conduct fit
        if length(range) == 1
            #Fit is always exact
            ys[i] = y[i]
        else
            p0::Array{Float64,1} = zeros(Float64,minimum([deg+1,length(range)-1]))
            fit = LsqFit.curve_fit(f,range-i,y[range],w,p0)
            ys[i] = fit.param[1]
        end
    end

    if !robust; return ys; end

    for iter=1:5
        #Find outliers
        resid::Array{Float64,1} = ys-y
        madev6::Float64 = 6*median(abs(resid))
        inliers::Array{Bool,1} = (abs(resid) .< madev6)

        #Compute Robust weights
        rw::Array{Float64,1} = zeros(Float64,length(y))
        rw[inliers] = (1-(resid/madev6).^2).^2

        #Conduct regression
        for i=1:length(y)
            #Get regression region
            if s > i-1
                range::Array{Int64,1} = [1:N]
            elseif s > length(y)-i
                range = [(length(y)-N+1):length(y)]
            else
                range = [(i-s):(i+s)]
            end

            #Compute local weights
            w = (1-(abs(range-i)/maximum(abs(range-i))).^3).^3

            #Conduct fit
            if length(range) == 1
                #Fit is always exact
                ys[i] = y[i]
            else
                p0::Array{Float64,1} = zeros(Float64,minimum([deg+1,length(range)-1]))
                fit = LsqFit.curve_fit(f,range-i,y[range],w.*rw[range],p0)
                ys[i] = fit.param[1]
            end
        end
    end

    return ys
end

ext: jl


text:
@testset "Model scoring and likelihood functions" begin
    #test a trivial scoring example
    obs=[BioSequences.LongSequence{DNAAlphabet{2}}("AAAAA")]
    order_seqs = BioBackgroundModels.get_order_n_seqs(obs, 0)
    coded_seqs = BioBackgroundModels.code_seqs(order_seqs)
    obs=Array(transpose(coded_seqs))

    source_pwm = zeros(1,4)
    source_pwm[1,1] = 1
    log_pwm = log.(source_pwm)

    source_stop=5

    #score prellocates
    ss_srcs = [revcomp_pwm(log_pwm)]
    ds_srcs=[cat(log_pwm,revcomp_pwm(log_pwm),dims=3)]
    score_mat_ds=zeros(source_stop,2)
    score_mat_ss=zeros(source_stop,1)
    score_matrices_ds=Vector{Matrix{Float64}}(undef, 2)
    score_matrices_ss=Vector{Vector{Float64}}(undef, 2)

    score_sources_ds!(score_mat_ds, score_matrices_ds, ds_srcs, obs[:,1], source_stop)
    @test score_matrices_ds[1] == [0.0 -Inf; 0.0 -Inf; 0.0 -Inf; 0.0 -Inf; 0.0 -Inf]

    score_sources_ss!(score_mat_ss, score_matrices_ss, ss_srcs, obs[:,1], source_stop)

    @test score_matrices_ss[1] == [-Inf for i in 1:5]

    #test a more complicated scoring example with two sources
    obs=[BioSequences.LongSequence{DNAAlphabet{2}}("ATGATGATGATG")]
    order_seqs = BioBackgroundModels.get_order_n_seqs(obs, 0)
    coded_seqs = BioBackgroundModels.code_seqs(order_seqs)
    obs=Array(transpose(coded_seqs))

    source_pwm = [.7 .1 .1 .1
                  .1 .1 .1 .7
                  .1 .1 .7 .1]

    source_pwm_2 = [.6 .1 .1 .2
    .2 .1 .1 .6
    .1 .2 .6 .1]

    target_s1 = [.7^3 .1^3; .1^3 .1^3; .1^3 (.1*.7^2); .7^3 .1^3; .1^3 .1^3; .1^3 (.1*.7^2); .7^3 .1^3; .1^3 .1^3; .1^3 (.1*.7^2); .7^3 .1^3]
    target_s2 = [.6^3 (.2*.1^2); (.2*.1^2) (.2*.1^2); (.2*.1^2) (.2*.6^2); .6^3 (.2*.1^2); (.2*.1^2) (.2*.1^2); (.2*.1^2) (.2*.6^2); .6^3 (.2*.1^2); (.2*.1^2) (.2*.1^2); (.2*.1^2) (.2*.6^2); .6^3 (.2*.1^2)]
    log_pwms = [log.(source_pwm), log.(source_pwm_2)]

    source_start = 1
    source_stop = 10

    #score prellocates
    ss_srcs = log_pwms
    ds_srcs=[cat(pwm,revcomp_pwm(pwm),dims=3) for pwm in log_pwms]
    score_mat_ds=zeros(source_stop,2)
    score_mat_ss=zeros(source_stop,1)

    score_sources_ds!(score_mat_ds, score_matrices_ds, ds_srcs, obs[:,1], [source_stop, source_stop])
    @test isapprox(score_matrices_ds[1], log.(target_s1))
    @test isapprox(score_matrices_ds[2], log.(target_s2))

    score_sources_ss!(score_mat_ss, score_matrices_ss, ss_srcs, obs[:,1], [source_stop, source_stop])
    @test isapprox(score_matrices_ss[1], log.(target_s1[:,1]))
    @test isapprox(score_matrices_ss[2], log.(target_s2[:,1]))

    #test score weaving and IPM likelihood calculations
    #trivial example using fake cardinality_penalty
    target=logaddexp(log(.25), 0.)
    lh_vec=zeros(2)
    obsl=1
    bg_scores=view([log(.25)],:,:)
    score_mat=[[0.]]
    osi=[1]
    source_wmls=[1]
    lme=0.
    cardinality_penalty=0. #not correct but for test purposes
    osi_emit=Vector{Int64}()

    @test weave_scores_ss!(lh_vec, obsl, bg_scores, score_mat, osi, source_wmls, lme, cardinality_penalty, osi_emit) == target

    target2=logsumexp([log(.25), log(.5), log(.5)])
    lh_vec=zeros(2)
    weavevec=zeros(3)
    score_mat=[zeros(1,2)]
    empty!(osi_emit)
    lme=log(.5)

    @test weave_scores_ds!(weavevec, lh_vec, obsl, bg_scores, score_mat, osi, source_wmls, lme, cardinality_penalty, osi_emit) == target2

    #test more complicated weaves
    obs=[BioSequences.LongSequence{DNAAlphabet{2}}("ATGATGATGATG")
    BioSequences.LongSequence{DNAAlphabet{2}}("TGATGATGATGA")]
    order_seqs = BioBackgroundModels.get_order_n_seqs(obs, 0)
    coded_seqs = BioBackgroundModels.code_seqs(order_seqs)
    obs=Array(transpose(coded_seqs))

    sources = [(log.(source_pwm), 1),(log.(source_pwm_2), 1)]
    dbl1_sources=[(log.(source_pwm), 1),(log.(source_pwm), 1)]

    ds_srcs=[cat(pwm[1],revcomp_pwm(pwm[1]),dims=3) for pwm in sources]
    dbl1_srcs=[cat(pwm[1],revcomp_pwm(pwm[1]),dims=3) for pwm in dbl1_sources]

    source_stops=[10,10]

    target_o2_s1 = [.1^3 .1^3; .1^3 (.1*.7^2); .7^3 .1^3; .1^3 .1^3; .1^3 (.1*.7^2); .7^3 .1^3; .1^3 .1^3; .1^3 (.1*.7^2); .7^3 .1^3; .1^3 .1^3]
    target_o2_s2 = [(.2*.1^2) (.2*.1^2); (.2*.1^2) (.2*.6^2); .6^3 (.2*.1^2); (.2*.1^2) (.2*.1^2); (.2*.1^2) (.2*.6^2); .6^3 (.2*.1^2); (.2*.1^2) (.2*.1^2); (.2*.1^2) (.2*.6^2); .6^3 (.2*.1^2); (.2*.1^2) (.2*.1^2)]

    score_mat_ds=zeros(source_stop,2)
    score_matrices_ds=Vector{Matrix{Float64}}(undef, 2)

    score_sources_ds!(score_mat_ds, score_matrices_ds, ds_srcs, obs[:,1], source_stops)

    o1_mats = copy(score_matrices_ds)

    @test isapprox(o1_mats, [log.(target_s1), log.(target_s2)])

    score_sources_ds!(score_mat_ds, score_matrices_ds, ds_srcs, obs[:,2], source_stops)

    o2_mats = copy(score_matrices_ds)

    @test isapprox(o2_mats, [log.(target_o2_s1), log.(target_o2_s2)])

    bg_scores = log.(fill(.5, (12,2)))
    log_motif_expectation = log(0.5 / size(bg_scores,1))
    osi = [1,2]
    obs_cardinality = length(osi)
    penalty_sum = sum(exp.(fill(log_motif_expectation,obs_cardinality)))
    cardinality_penalty=log(1.0-penalty_sum)
    empty!(osi_emit)
    source_wmls=[3,3]

    lh_target = -8.87035766177774
    lh_vec=zeros(13)

    o1_lh = weave_scores_ds!(weavevec, lh_vec, 12, view(bg_scores,:,1), o1_mats, osi, source_wmls, log_motif_expectation, cardinality_penalty, osi_emit)
    @test isapprox(lh_target,o1_lh)

    empty!(osi_emit)

    o2_lh = weave_scores_ds!(weavevec, lh_vec, 12, view(bg_scores,:,2), o2_mats, osi, source_wmls, log_motif_expectation, cardinality_penalty, osi_emit)

    lh,cache = IPM_likelihood(sources, obs, [12,12], bg_scores, trues(2,2), true,true)
    @test isapprox(o1_lh,cache[1])
    @test isapprox(o2_lh,cache[2])
    @test isapprox(lps(o1_lh,o2_lh),lh)

    #test source penalization
    score_sources_ds!(score_mat_ds, score_matrices_ds, dbl1_srcs, obs[:,1], source_stops)

    dbl_score_mat=copy(score_matrices_ds)

    empty!(osi_emit)

    naive=weave_scores_ds!(weavevec, lh_vec, 12, view(bg_scores,:,1),Vector{Matrix{Float64}}(), Vector{Int64}(), Vector{Int64}(), log_motif_expectation, cardinality_penalty, osi_emit)

    empty!(osi_emit)

    single=weave_scores_ds!(weavevec, lh_vec, 12, view(bg_scores,:,1), [dbl_score_mat[1]], [1], [3], log_motif_expectation, cardinality_penalty, osi_emit)

    empty!(osi_emit)

    double=weave_scores_ds!(weavevec, lh_vec, 12, view(bg_scores,:,1), dbl_score_mat, [1,2], [3,3], log_motif_expectation, cardinality_penalty, osi_emit)

    empty!(osi_emit)

    triple=weave_scores_ds!(weavevec, lh_vec, 12, view(bg_scores,:,1), [dbl_score_mat[1] for i in 1:3], [1,2,3], [3,3,3], log_motif_expectation, cardinality_penalty, osi_emit)

    @test (single-naive) > (double-single) > (triple-double)

    naive_target=-16.635532333438686
    naive = IPM_likelihood(sources, obs, [findfirst(iszero,obs[:,o])-1 for o in 1:size(obs)[2]], bg_scores,falses(size(obs)[2],length(sources)))
    @test naive==naive_target

    #test IPM_likelihood clean vector and cache calculations
    mix_matrix=trues(2,2)
    baselh,basecache = IPM_likelihood(sources, obs, [12,12], bg_scores, mix_matrix, true,true)

    clean=[true,false]

    unchangedlh,unchangedcache=IPM_likelihood(sources, obs, [12,12], bg_scores, mix_matrix, true, true, basecache, clean)

    changed_lh,changedcache=IPM_likelihood(sources, obs, [12,12], bg_scores, BitMatrix([true true; false false]), true, true, basecache, clean)

    indep_lh, indepcache=IPM_likelihood(sources, obs,[12,12], bg_scores, BitMatrix([true true; false false]), true, true)

    @test baselh==unchangedlh!=changed_lh==indep_lh
    @test basecache==unchangedcache!=changedcache==indepcache

    #check that IPM_likelihood works in single stranded operation
    IPM_likelihood(sources, obs,[12,12], bg_scores, BitMatrix([true true; false false]), false)
end


text:
# run last
using Bridge, Test, StaticArrays

Random.seed!(1)
tt = collect(range(0.0, stop=1.0, length=50))
W = sample(tt, Wiener{SVector{3,Float64}}())
VW = VSamplePath(tt, copy(Bridge.mat(W.yy)))
Random.seed!(1)
sample!(VW, Wiener())

@test maximum(abs.(VW.yy - Bridge.mat(W.yy))) < eps()


text:
module CardiovascularModels

using ModelingToolkit, OrdinaryDiffEq

@parameters t
D = Differential(t)

include("abstractComponents.jl")
export Con, CompliantCompartment, InertiaCompartment, VolumelessComponent, PressurizedCompartment
include("physiologicalComponents.jl")
export Ventricle, ElasticVessel, DampedElasticVessel
export qrsDriver, automaticDriver, ConstantPressure, Resistor, Valve

include("helpers.jl")

export serial_connect

end # module


text:
if !@isdefined(TabLevel)
    TabLevel = ""
end
println(TabLevel*"Open MajoranaWolfReps.jl")
TabLevel=TabLevel*"    "

if HasMathLink ###Extra Math Link dependents

    using MathLinkExtras

    function ScaleToStr(scale::Union{MathLink.WExpr,MathLink.WSymbol},tex::Bool,KetState::Bool,NumOps::Int,FirstTerm::Bool)
        if FirstTerm
            return "$(scale)*"
        else
            return "+ $(scale)*"
        end
    end

end

TabLevel=TabLevel[1:end-4]
println(TabLevel*"Close MajoranaWolfReps.jl")


text:
#!/usr/bin/env julia

include("goc_challenge1_huristic.jl")
include("scopf-main.jl")

scopf_main(parse_scopf_commandline())


text:
# This file contains all functions that either do not directly interface with `brew`
# or peek into the internals of Homebrew a bit, such as `install_brew()` which
# installs Homebrew using Git, `tap_exists()` which presupposes knowledge of the
# internal layout of Taps, and `installed()` which presupposed knowledge of the
# location of files within the Cellar

# This is the tap that contains our manually-curated overrides
const tapname = "staticfloat/juliadeps"
const tappath = joinpath(brew_prefix,"Library","Taps","staticfloat","homebrew-juliadeps")

# This is the tap that will contain our automatically-translated overrides
const auto_tapname = "staticfloat/juliatranslated"
const auto_tappath = joinpath(brew_prefix,"Library","Taps","staticfloat","homebrew-juliatranslated")

# Where we download brew from
const BREW_URL = "https://github.com/Homebrew/brew"
const BREW_BRANCH = "master"
const BREW_STABLE_SHA = "85e4013989e849b374ff40858969be1125d63370"
"""
`install_brew()`

Ensures that Homebrew is installed as desired, that our basic Taps are available
and that we have whatever binary tools we need, such as `install_name_tool`
"""
function install_brew()
    # Ensure brew_prefix exists
    if !isdir(brew_prefix) || !isfile(joinpath(brew_prefix,"bin","brew"))
        try
            mkdir(brew_prefix)
        catch
        end

        try
            @info("Downloading brew...")
            run(pipeline(`curl -\# -L $BREW_URL/tarball/$BREW_BRANCH`,
                                 `tar xz -m --strip 1 -C $brew_prefix`))
        catch
            @warn("Could not download/extract $BREW_URL/tarball/$BREW_BRANCH into $(brew_prefix)!")
            rethrow()
        end

        cd(brew_prefix) do
            try
                git(`fetch --unshallow`)
            catch
            end
        end
    end

    # Tap homebrew/core, always and forever
    tap("homebrew/core")

    # Tap our own "overrides" taps
    tap("staticfloat/juliadeps")
    tap("staticfloat/juliatranslated")

    if !clt_installed() && !installed("cctools")
        # If we don't have the command-line tools installed, then let's grab
        # cctools, as we need that to install bottles that need relocation
        add("cctools")
    end

    if !git_installed() && !installed("git")
        # If we don't have a git available, install it now
        add("git")
    end

    # Update the environment before running update()
    update_env()

    # Update immediately so that we get onto our "stable" tag
    update(verbose=true)

    return nothing
end

"""
`tap_exists(tap_name::AbstractString)`

Check to see if a tap called `tap_name` (ex: `"staticfloat/juliadeps"`) exists
"""
function tap_exists(tap_name::AbstractString)
    path = joinpath(brew_prefix,"Library","Taps", dirname(tap_name), "homebrew-$(basename(tap_name))")
    return isdir(path)
end

"""
`installed(pkg::Union{AbstractString,BrewPkg})`

Return true if the given package `pkg` is a directory in the Cellar, showing
that it has been installed (but possibly not linked, see `linked()`)
"""
function installed(pkg::StringOrPkg) end

function installed(name::AbstractString)
    isdir(joinpath(brew_prefix,"Cellar",basename(name)))
end

function installed(pkg::BrewPkg)
    installed(pkg.name)
end

"""
`linked(pkg::Union{AbstractString,BrewPkg})`

Returns true if the given package `pkg` is linked to LinkedKegs, signifying
all files installed by this package have been linked into the global prefix.
"""
function linked(pkg::StringOrPkg) end

function linked(name::AbstractString)
    return islink(joinpath(brew_prefix,"var","homebrew","linked",basename(name)))
end

function linked(pkg::BrewPkg)
    return linked(pkg.name)
end

"""
`formula_path(pkg::Union{AbstractString,BrewPkg})`

Returns the absolute path on-disk of the given package `pkg`.
"""
function formula_path(pkg::StringOrPkg) end

function formula_path(name::AbstractString)
    path, tap_path = formula_tap(name)

    if isempty(tap_path)
        return joinpath(brew_prefix, "Library", "Taps", "homebrew", "homebrew-core", "Formula", "$path.rb")
    else
        # Insert the "homebrew-" that exists in all taps
        tap_path = "$(dirname(tap_path))/homebrew-$(basename(tap_path))"
        return joinpath(brew_prefix, "Library", "Taps", tap_path, "$path.rb")
    end
end

function formula_path(pkg::BrewPkg)
    return formula_path(fullname(pkg))
end

"""
`update_tag(;verbose::Bool = false)`

We maintain our own "stable" tag that overrides Homebrew so that we can
update at our own pace along with them.  Make sure to call `update_env()`
before calling this, as calling our `git` doesn't work properly otherwise.
"""
function update_tag(;verbose::Bool=false)
    global BREW_STABLE_SHA, brew_prefix

    # If a recent enough version of git is not installed, then install our own
    # This is necessary because some people do not have git installed but already
    # have Homebrew installed, and the check up above during install_brew()
    # doesn't run every time
    if !git_installed()
        add("git")
    end

    if verbose
        git = cmd -> run(`git $cmd`)
    else
        git = cmd -> run(pipeline(`git $cmd`, stdout=devnull, stderr=devnull))
    end

    cd(brew_prefix) do
        try
            git(`tag -d 9.9.9`)
        catch
        end
        git(`tag 9.9.9 $BREW_STABLE_SHA`)
        git(`checkout 9.9.9`)
    end
    return nothing
end


text:
#!/usr/bin/env julia

# Dependencies
# 7zip (Windows)

include("jl/julia-download.jl")

# Julia Packages

pkgname = "pkg"^20
mkpath("julia/pkg")
mv("julia/pkg", "julia/$pkgname")

try
  ENV["JULIA_PKGDIR"] = abspath("julia/$pkgname")
  @osx_only run(`julia/julia-$ver-mac/bin/julia -i --color=yes jl/packages.jl`)
  @windows_only run(`julia/julia-$ver-win64/bin/julia.exe -i jl/packages.jl`)
  delete!(ENV, "JULIA_PKGDIR")
finally
  mv("julia/$pkgname", "julia/pkg")
end

# Atom Build

isdir("out/juno") && rm("out/juno", recursive=true)
isdir("out/installer") && rm("out/installer", recursive=true)

run(`node script/build --build-dir=out --no-install --channel=stable
                       --no-auto-update --julia-version=$ver`)

max_stars_repo_name: JunoLab/atom

max_issues_repo_name: JunoLab/atom

max_forks_repo_name: JunoLab/atom


text:
module JOLI4Flux

    using JOLI: joAbstractLinearOperator
    using ChainRulesCore, FillArrays
    import Base.*

    function ChainRulesCore.rrule(::typeof(*), A::joAbstractLinearOperator{ADDT,ARDT}, v) where {ADDT,ARDT}
      Y = A*v
      function time_pullback(dy)
        DY = unthunk(dy)
        return NoTangent(), NoTangent(), @thunk(A' * DY)
      end
      return Y, time_pullback
    end

    function *(A::joAbstractLinearOperator{ADDT,ARDT}, v::FillArrays.Fill{vDT,1}) where {ADDT,ARDT,vDT<:Number}
        A*Vector(v)
    end

end


text:
"""
    AbstractWord{T} <: AbstractVector{T}
Abstract type representing words over an Alphabet.

`AbstractWord` as such has its meaning only in the contex of an Alphabet.
The subtypes of `AbstractWord{T}` need to implement the following methods which
constitute `AbstractWord` interface:
 * a constructor from `AbstractVector{T}`
 * linear indexing (1-based) consistent with iteration returning pointers to letters of an alphabet (`getindex`, `setindex`, `length`),
 * `length`: the length of word as written in the alphabet,
 * `Base.push!`/`Base.pushfirst!`: append a single value at the end/beginning,
 * `Base.pop!`/`Base.popfirst!`: pop a single value from the end/beginning,
 * `Base.append!`/`Base.prepend!`: append a another word at the end/beginning,
 * `Base.resize!`: drop/extend a word at the end to the requested length
 * `Base.:*`: word concatenation (monoid binary operation),
 * `Base.similar`: an uninitialized word of a similar type/storage.

Note that `length` represents free reduced word (how it is written in an alphabet)
and not its the shortest form (e.g. the normal form).

The following are implemented for `AbstractWords` but can be overloaded for
performance reasons:

* `Base.==`: the equality (as words),
* `Base.hash`: simple uniqueness hashing function
* `Base.view`: creating `SubWord` e.g. based on subarray.
"""
abstract type AbstractWord{T<:Integer} <: AbstractVector{T} end

function Base.hash(w::AbstractWord, h::UInt)
    h = hash(AbstractWord, h)
    for i in w
        h = hash(i, h)
    end
    return h
    # foldl((h, x) -> hash(x, h), w, init = hash(AbstractWord, h))
end

@inline Base.:(==)(w::AbstractWord, v::AbstractWord) =
    length(w) == length(v) && all(@inbounds w[i] == v[i] for i in 1:length(w))

Base.convert(::Type{W}, w::AbstractWord) where W<:AbstractWord = W(w)
Base.convert(::Type{W}, w::W) where W<:AbstractWord = w

# resize + copyto!
function store!(w::AbstractWord, v::AbstractWord)
    resize!(w, length(v))
    copyto!(w, v)
    return w
end

Base.size(w::AbstractWord) = (length(w),)

Base.one(::Type{W}) where {T, W<:AbstractWord{T}} = W(T[])
Base.one(::W) where W <: AbstractWord = one(W)
Base.isone(w::AbstractWord) = iszero(length(w))

Base.getindex(w::W, u::AbstractRange) where W<:AbstractWord =
    W([w[i] for i in u])

Base.view(w::AbstractWord, u::AbstractRange) = w[u] # general fallback

Base.:^(w::AbstractWord, n::Integer) = n >= 0 ? Base.power_by_squaring(w, n) :
    throw(DomainError(n, "To rise a Word to negative power you need to provide its inverse."))
Base.literal_pow(::typeof(^), w::AbstractWord, ::Val{p}) where p =
    p >= 0 ? Base.power_by_squaring(w, p) :
    throw(DomainError(p, "To rise a Word to negative power you need to provide its inverse."))

function Base.findnext(subword::AbstractWord, word::AbstractWord, pos::Integer)
    k = length(subword)
    f = first(subword)
    @inbounds for i in pos:length(word)-k+1
        word[i] == f || continue
        issub = true
        for j in 2:k
            if word[i+j-1] != subword[j]
                issub = false
                break
            end
        end
        issub == true && return i:i+k-1
    end
    return nothing
end

@inline Base.findfirst(subword::AbstractWord, word::AbstractWord) = findnext(subword, word, firstindex(word))
@inline Base.occursin(subword::AbstractWord, word::AbstractWord) = findfirst(subword, word) !== nothing

"""
    longestcommonprefix(u::AbstractWord, v::AbstractWord)
Returns the length of longest common prefix of two words (and simultaneously
the index at which the prefix ends).
"""
function longestcommonprefix(u::AbstractWord, v::AbstractWord)
    n=0
    for (lu, lv) in zip(u,v)
        lu != lv && break
        n += 1
    end
    return n
end

"""
    lcp(u::AbstractWord, v::AbstractWord)
See [`longestcommonprefix`](@ref).
"""
lcp(u::AbstractWord, v::AbstractWord) = longestcommonprefix(u,v)

"""
    isprefix(u::AbstractWord, v::AbstractWord)
Check if `u` is a prefix of `v`.
"""
@inline function isprefix(u::AbstractWord, v::AbstractWord)
    k = length(u)
    k < length(v) || return false
    @inbounds for i in eachindex(u)
        u[i] == v[i] || return false
    end
    return true
end

"""
    issuffix(u::AbstractWord, v::AbstractWord)
Check if `u` is a suffix of `v`.
"""
@inline function issuffix(u::AbstractWord, v::AbstractWord)
    k = length(u)
    k ≤ length(v) || return false
    @inbounds for i in eachindex(u)
        u[i] == v[end-k+i] || return false
    end
    return true
end

function Base.show(io::IO, ::MIME"text/plain", w::AbstractWord)
    print(io, typeof(w), ": ")
    show(io, w)
end

function Base.show(io::IO, w::AbstractWord{T}) where T
    if isone(w)
        print(io, "(id)")
    else
        join(io, w, "·")
    end
end

_max_alphabet_length(::Type{<:AbstractWord{T}}) where T = typemax(T)


text:
module FefferyAntdComponents
using Dash

const resources_path = realpath(joinpath( @__DIR__, "..", "deps"))
const version = "0.0.1rc3"

include("jl/'feffery'_antdalert.jl")
include("jl/'feffery'_antdanchor.jl")
include("jl/'feffery'_antdbutton.jl")
include("jl/'feffery'_antdcascader.jl")
include("jl/'feffery'_antdcheckbox.jl")
include("jl/'feffery'_antdcheckboxgroup.jl")
include("jl/'feffery'_antdcollapse.jl")
include("jl/'feffery'_antddatepicker.jl")
include("jl/'feffery'_antddaterangepicker.jl")
include("jl/'feffery'_antddivider.jl")
include("jl/'feffery'_antddrawer.jl")
include("jl/'feffery'_antdempty.jl")
include("jl/'feffery'_antdicon.jl")
include("jl/'feffery'_antdinput.jl")
include("jl/'feffery'_antdmenu.jl")
include("jl/'feffery'_antdmessage.jl")
include("jl/'feffery'_antdmodal.jl")
include("jl/'feffery'_antdnotification.jl")
include("jl/'feffery'_antdpagination.jl")
include("jl/'feffery'_antdpopover.jl")
include("jl/'feffery'_antdradiogroup.jl")
include("jl/'feffery'_antdresult.jl")
include("jl/'feffery'_antdselect.jl")
include("jl/'feffery'_antdskeleton.jl")
include("jl/'feffery'_antdslider.jl")
include("jl/'feffery'_antdspace.jl")
include("jl/'feffery'_antdspin.jl")
include("jl/'feffery'_antdsteps.jl")
include("jl/'feffery'_antdswitch.jl")
include("jl/'feffery'_antdtable.jl")
include("jl/'feffery'_antdtag.jl")
include("jl/'feffery'_antdtooltip.jl")
include("jl/'feffery'_antdtransfer.jl")
include("jl/'feffery'_antdtree.jl")
include("jl/'feffery'_antdtreeselect.jl")
include("jl/'feffery'_antdupload.jl")
include("jl/'feffery'_antdcol.jl")
include("jl/'feffery'_antdrow.jl")
include("jl/'feffery'_antdcontent.jl")
include("jl/'feffery'_antdfooter.jl")
include("jl/'feffery'_antdheader.jl")
include("jl/'feffery'_antdlayout.jl")
include("jl/'feffery'_antdsider.jl")
include("jl/'feffery'_antdtabpane.jl")
include("jl/'feffery'_antdtabs.jl")
include("jl/'feffery'_antdparagraph.jl")
include("jl/'feffery'_antdtext.jl")
include("jl/'feffery'_antdtitle.jl")

function __init__()
    DashBase.register_package(
        DashBase.ResourcePkg(
            "feffery_antd_components",
            resources_path,
            version = version,
            [
                DashBase.Resource(
    relative_package_path = "feffery_antd_components.min.js",
    external_url = "https://unpkg.com/feffery_antd_components@0.0.1rc3/feffery_antd_components/feffery_antd_components.min.js",
    dynamic = nothing,
    async = nothing,
    type = :js
),
DashBase.Resource(
    relative_package_path = "feffery_antd_components.min.js.map",
    external_url = "https://unpkg.com/feffery_antd_components@0.0.1rc3/feffery_antd_components/feffery_antd_components.min.js.map",
    dynamic = true,
    async = nothing,
    type = :js
)
            ]
        )

    )
end
end

max_stars_repo_name: Beddow/feffery-antd-components

max_issues_repo_name: Beddow/feffery-antd-components

max_forks_repo_name: Beddow/feffery-antd-components


text:
"""
Original Kings law (1910) represents the relation between voltage and velocity
in a hotwire anemometer.
The law is given by:

`E^2 = A + B*U^0.5`

This function estimates `A` and `B`.
"""
linear_king_fit(E, U) = linear_fit(sqrt.(U), E .* E)

"Type that represents a Linear (original) King's law"
struct LinearKingFit{T<:Number} <: AbstractLeastSquares
    coefs::NTuple{2,T} #Array{Float64,1}
    LinearKingFit(A::T, B::T) where {T<:Number} = new{T}((A,B))
    LinearKingFit(coefs::NTuple{2,T}) where {T<:Number}= new{T}(coefs)
    LinearKingFit(E::AbstractVector{T},U::AbstractVector{T}) where {T <: Number} = new{T}(linear_king_fit(E, U))
end

(f::LinearKingFit)(E) = ( (E.*E .- f.coefs[1]) ./ f.coefs[2] ) .^ 2

"Equation that computes the error of the modified King's law "
kingfun(x, a) = a[1] + a[2] * x[2] ^ a[3] - x[1]*x[1]

"""
Uses nonlinear least squares to fit the modified King's law:

`E^2 = A + B*U^n`

The Original (linear) King's law is used to estimate `A` and `B` when `n=1/2`.
This initial value is used as an initial guess for fitting the nonlinear modified King's law
using the function `nonlinear_fit`.
"""
function king_fit(E, U, eps=1e-8, maxiter=200)

    f = linear_king_fit(E, U)
    a0 = [f[1], f[2], 0.5]

    coefs, converged, niter = nonlinear_fit(hcat(E,U), kingfun, a0, eps, maxiter)
    return (coefs[1], coefs[2], coefs[3])
end

"Type that represents the modified King's law "
struct KingFit{T<:Number} <: AbstractLeastSquares
    coefs::NTuple{3,T} #Array{Float64,1}
    KingFit(A::T, B::T, n::T) where {T<:Number} = new{T}((A, B, n)) #copy(coefs))
    KingFit(coefs::NTuple{3,T}) where {T<:Number} = new{T}(coefs)
    KingFit(A::T, B::T) where {T<:Number} = new{T}((A, B, one(T)/2))
    KingFit(coefs::NTuple{2,T}) where {T<:Number} = new{T}(coefs[1], coefs[2])
end
KingFit(E,U, eps=1e-8, maxiter=200) = KingFit(king_fit(E, U, eps, maxiter))

(f::KingFit)(E) = ( (E*E - f.coefs[1]) / f.coefs[2]) ^ (1 / f.coefs[3])


text:
### math.jl
#
# Copyright (C) 2016, 2017 Mosè Giordano.
#
# Maintainer: Mosè Giordano <mose AT gnu DOT org>
# Keywords: uncertainty, error propagation, physics
#
# This file is a part of Measurements.jl.
#
# License is MIT "Expat".
#
### Commentary:
#
# This file contains definition of mathematical functions that support
# Measurement objects.
#
# Note: some functions defined here (like all degree-related and reciprocal
# trigonometric functions, fld, cld, hypot, cbrt, abs, mod) are redundant in the
# sense that you would get the correct result also without their definitions,
# but having them defined here avoids some calculations and slightly improves
# performance.  Likewise, multiple methods are provided for functions taking two
# (or more) arguments because when only one argument is of Measurement type we
# can use the simple `result' function for one derivative that is faster than
# the generic method.
#
### Code:

export @uncertain

# This function is to be used by methods of mathematical operations to produce a
# `Measurement' object in output.  Arguments are:
#   * val: the nominal result of operation G(a)
#   * der: the derivative ∂G/∂a of G with respect to the variable a
#   * a: the only argument of G
# In this simple case of unary function, we don't have the problem of correlated
# variables (thus making this method much faster than the next one), so we can
# calculate the uncertainty of G(a) as
#   σ_G = |σ_a·∂G/∂a|
# The list of derivatives with respect to each measurement is updated with
#   ∂G/∂a · previous_derivatives
@inline function result(val::T, der::Real, a::Measurement{<:AbstractFloat}) where {T<:Real}
    newder = empty_der1(a)
    @inbounds for tag in keys(a.der)
        if ! iszero(tag[2]) # Skip values with 0 uncertainty
            newder = Derivatives(newder, tag=>der*a.der[tag])
        end
    end
    # If uncertainty of "a" is null, the uncertainty of result is null as well,
    # even if the derivative is NaN or infinite.  In any other case, use
    # σ_G = |σ_a·∂G/∂a|.
    σ = iszero(a.err) ? a.err : abs(der*a.err)
    # The tag for derived quantities is 0, for independent ones tag > 0.
    Measurement{float(T)}(val, σ, UInt64(0), newder)
end

# Get the common type parameter of a collection of Measurement objects.  The first two
# methods are for the trivial cases of homogeneous tuples and arrays, the last, inefficient,
# method is for inhomogeneous collections (probably the least common case).
gettype(::Tuple{Vararg{Measurement{T}}}) where {T<:AbstractFloat} = T
gettype(::AbstractArray{Measurement{T}}) where {T<:AbstractFloat} = T
_eltype(::Measurement{T}) where {T<:AbstractFloat} = T
gettype(collection) = promote_type(_eltype.(collection)...)

# This function is similar to the previous one, but applies to mathematical
# operations with more than one argument, so the formula to propagate
# uncertainty is more complicated because we have to take into account
# correlation between arguments.  The arguments are the same as above, but `der'
# and `a' are tuples of the same length (`der' has the derivatives of G with
# respect to the corresponding variable in `a').
#
# Suppose we have a function G = G(a1, a2) of two arguments.  a1 and a2 are
# correlated, because they come from some mathematical operations on really
# independent variables x, y, z, say a1 = a1(x, y), a2 = a2(x, z).  The
# uncertainty on G(a1, a2) is calculated as follows:
#   σ_G = sqrt((σ_x·∂G/∂x)^2 + (σ_y·∂G/∂y)^2 + (σ_z·∂G/∂z)^2)
# where ∂G/∂x is the partial derivative of G with respect to x, and so on.  We
# can expand the previous formula to:
#   σ_G = sqrt((σ_x·(∂G/∂a1·∂a1/∂x + ∂G/∂a2·∂a2/∂x))^2 + (σ_y·∂G/∂a1·∂a1/∂y)^2 +
#               + (σ_z·∂G/∂a2·∂a2/∂z)^2)
@inline function result(val, der, a)
    @assert length(der) == length(a)
    T = gettype(a)
    nil::T = zero(T)
    err = nil
    newder = empty_der2(nil)
    # Iterate over all independent variables.  We first iterate over all
    # variables listed in `a' in order to get all independent variables upon
    # which those variables depend, then we get the `tag' of each independent
    # variable, skipping variables that have been already taken into account.
    @inbounds for y in a
        for tag in keys(y.der)
            if tag ∉ keys(newder) # Skip independent variables already considered
                σ_x = tag[2]
                if ! iszero(σ_x)  # Skip values with 0 uncertainty
                    ∂G_∂x::T = nil
                    # Iteratate over all the arguments of the function
                    for (i, x) in enumerate(a)
                        # Calculate the derivative of G with respect to the
                        # current independent variable.  In the case of the x
                        # independent variable of the example above, we should
                        # get   ∂G/∂x = ∂G/∂a1·∂a1/∂x + ∂G/∂a2·∂a2/∂x
                        ∂a_∂x = derivative(x, tag) # ∂a_i/∂x
                        if ! iszero(∂a_∂x) # Skip values with 0 partial derivative
                            # der[i] = ∂G/∂a_i
                            ∂G_∂x = ∂G_∂x + der[i]*∂a_∂x
                        end
                    end
                    if ! iszero(∂G_∂x)
                        # Add (σ_x·∂G/∂x)^2 to the total uncertainty (squared), but only if
                        # the derivative is not zero.
                        newder = Derivatives(newder, tag=>∂G_∂x)
                        err = err + abs2(σ_x*∂G_∂x)
                    end
                end
            end
        end
    end
    return Measurement(T(val), sqrt(err), UInt64(0), newder)
end

# "result" function for complex-valued functions of one real argument (like "besselh").
# This takes the same argument as the first implementation of "result", but with complex
# "val" and "der".
result(val::Complex, der::Complex, a::Measurement) =
    complex(result(real(val), real(der), a), result(imag(val), imag(der), a))

# "result" function for complex-valued functions of one complex argument:
#
#   f(z) = f(x, y) = u(x, y) + im * v(x, y)
#
# where z = x + im * y, x and y are real variables, u and v are real-valued functions.
# Arguments of this methods are:
#
#   1) the nominal value of the result (a complex number)
#   2) the 4-tuple of derivatives (∂u/∂x, ∂u/∂y, ∂v/∂x, ∂v/∂y)
#   3) the input measurement
result(val::Complex, der, a::Complex{<:Measurement}) =
    complex(result(real(val), der[1:2], reim(a)),
            result(imag(val), der[3:4], reim(a)))

### @uncertain macro.
"""
    @uncertain f(value ± stddev, ...)

A macro to calculate `f(value) ± uncertainty`, with `uncertainty` derived from `stddev`
according to rules of linear error propagation theory.

Function `f` can accept any number of real arguments.
"""
macro uncertain(expr::Expr)
    f = esc(expr.args[1]) # Function name
    n = length(expr.args) - 1
    if n == 1
        a = esc(expr.args[2]) # Argument, of Measurement type
        return quote
            let x = measurement($a)
                result($f(x.val), Calculus.derivative($f, x.val), x)
            end
        end
    else
        a = expr.args[2:end] # Arguments, as an array of expressions
        args = :([])  # Build up array of arguments
        [push!(args.args, :(measurement($(esc(a[i]))))) for i=1:n] # Fill the array
        argsval =:([])  # Build up the array of values of arguments
        [push!(argsval.args, :($(args.args[i]).val)) for i=1:n] # Fill the array
        return :( result($f($argsval...),
                         Calculus.gradient(x -> $f(x...), $argsval),
                         $args) )
    end
end

### Elementary arithmetic operations:

# Addition: +
Base.:+(a::Measurement, b::Measurement) = result(a.val + b.val, (1, 1), (a, b))
Base.:+(a::Real, b::Measurement) = result(a + b.val, 1, b)
Base.:+(a::Measurement, b::Bool) = result(a.val + b, 1, a)
Base.:+(a::Measurement, b::Real) = result(a.val + b, 1, a)

# Subtraction: -
Base.:-(a::Measurement) = result(-a.val, -1, a)
Base.:-(a::Measurement, b::Measurement) = result(a.val - b.val, (1, -1), (a, b))
Base.:-(a::Real, b::Measurement) = result(a - b.val, -1, b)
Base.:-(a::Measurement, b::Real) = result(a.val - b, 1, a)

# Multiplication: *
function Base.:*(a::Measurement, b::Measurement)
    aval = a.val
    bval = b.val
    return result(aval*bval, (bval, aval), (a, b))
end
Base.:*(a::Bool, b::Measurement) = result(a*b.val, a, b)
Base.:*(a::Real, b::Measurement) = result(a*b.val, a, b)
Base.:*(a::Measurement, b::Bool) = result(a.val*b, b, a)
Base.:*(a::Measurement, b::Real) = result(a.val*b, b, a)

# muladd and fma

for f in (:fma, :muladd)
    @eval begin
        # All three arguments are Measurement
        function Base.$f(a::Measurement, b::Measurement, c::Measurement)
            x = a.val
            y = b.val
            z = c.val
            return result(($f)(x, y, z), (y, x, one(z)), (a, b, c))
        end

        # First argument is always Measurement
        function Base.$f(a::Measurement, b::Measurement, c::Real)
            x = a.val
            y = b.val
            return result(($f)(x, y, c), (y, x), (a, b))
        end

        function Base.$f(a::Measurement, b::Real, c::Measurement)
            x = a.val
            z = c.val
            return result(($f)(x, b, z), (b, one(z)), (a, c))
        end

        Base.$f(a::Measurement, b::Real, c::Real) =
            result(($f)(a.val, b, c), b, a)

        # Secon argument is always Measurement
        function Base.$f(a::Real, b::Measurement, c::Measurement)
            y = b.val
            z = c.val
            return result(($f)(a, y, z), (a, one(z)), (b, c))
        end

        Base.$f(a::Real, b::Measurement, c::Real) =
            result(($f)(a, b.val, c), a, b)

        # Third argument is Measurement
        function Base.$f(a::Real, b::Real, c::Measurement)
            z = c.val
            return result(($f)(a, b, z), one(z), c)
        end
    end
end

# Division: /, div, fld, cld
function Base.:/(a::Measurement, b::Measurement)
    x = a.val
    y = b.val
    oneovery = 1 / y
    return result(x / y, (oneovery, -x * abs2(oneovery)), (a, b))
end
Base.:/(a::Real, b::Measurement) = result(a/b.val, -a/abs2(b.val), b)
Base.:/(a::Measurement{T}, b::Real) where {T<:AbstractFloat} = result(a.val/b, 1/T(b), a)

# 0.0 as partial derivative for both arguments of "div", "fld", "cld" should be
# correct for most cases.  This has been tested against "@uncertain" macro.
Base.div(a::Measurement, b::Measurement) = result(div(a.val, b.val), (0, 0), (a, b))
Base.div(a::Measurement, b::Real) = result(div(a.val, b), 0, a)
Base.div(a::Real, b::Measurement) = result(div(a, b.val), 0, b)

Base.fld(a::Measurement, b::Measurement) = result(fld(a.val, b.val), (0, 0), (a, b))
Base.fld(a::Measurement, b::Real) = result(fld(a.val, b), 0, a)
Base.fld(a::Real, b::Measurement) = result(fld(a, b.val), 0, b)

Base.cld(a::Measurement, b::Measurement) = result(cld(a.val, b.val), (0, 0), (a, b))
Base.cld(a::Measurement, b::Real) = result(cld(a.val, b), 0, a)
Base.cld(a::Real, b::Measurement) = result(cld(a, b.val), 0, b)

# Inverse: inv
function Base.inv(a::Measurement)
    inverse = inv(a.val)
    return result(inverse, -abs2(inverse), a)
end

# signbit

Base.signbit(a::Measurement) = signbit(a.val)

# Power: ^

function Base.:^(a::Measurement, b::Measurement)
    aval = a.val
    bval = b.val
    pow = aval^bval
    return result(pow, (aval^(bval - 1)*bval, pow*log(aval)), (a, b))
end

function Base.:^(a::Measurement, b::Integer)
    x = a.val
    return result(x ^ b, b * x ^ (b - 1), a)
end

function Base.:^(a::Measurement{T},  b::Rational) where {T<:AbstractFloat}
    x = a.val
    return result(x ^ b, b * x ^ (b - one(T)), a)
end

function Base.:^(a::Measurement,  b::Real)
    x = a.val
    return result(x ^ b, b * x ^ (b - 1), a)
end

Base.:^(::Irrational{:ℯ}, b::Measurement) = exp(b)

function Base.:^(a::Real,  b::Measurement)
    res = a^b.val
    return result(res, res*log(a), b)
end

function Base.exp2(a::Measurement{T}) where {T<:AbstractFloat}
    pow = exp2(a.val)
    return result(pow, pow*log(T(2)), a)
end

### Trigonometric functions

# deg2rad, rad2deg

Base.deg2rad(a::Measurement) = a * (oftype(a.val, pi) / 180)
Base.rad2deg(a::Measurement) = a * (180 / oftype(a.val, pi))

# Cosine: cos, cosd, cosh

function Base.cos(a::Measurement)
    s, c = sincos(a.val)
    return result(c, -s, a)
end

function Base.cosd(a::Measurement)
    aval = a.val
    return result(cosd(aval), -deg2rad(sind(aval)), a)
end

function Base.cosh(a::Measurement)
    aval = a.val
    result(cosh(aval), sinh(aval), a)
end

# Sine: sin, sind, sinh

function Base.sin(a::Measurement)
    s, c = sincos(a.val)
    return result(s, c, a)
end

function Base.sind(a::Measurement)
    aval = a.val
    return result(sind(aval), deg2rad(cosd(aval)), a)
end

function Base.sinh(a::Measurement)
    aval = a.val
    result(sinh(aval), cosh(aval), a)
end

# Sincos: sincos

function Base.sincos(a::Measurement)
    s, c = sincos(a.val)
    return (result(s, c, a), result(c, -s, a))
end

if isdefined(Base, :sincospi)
    function Base.sincospi(a::Measurement)
        s, c = sincospi(a.val)
        return (result(s, c * pi, a), result(c, -s * pi, a))
    end
end

# Tangent: tan, tand, tanh

function Base.tan(a::Measurement)
    aval = a.val
    return result(tan(aval), abs2(sec(aval)), a)
end

function Base.tand(a::Measurement)
    aval = a.val
    return result(tand(aval), deg2rad(abs2(secd(aval))), a)
end

function Base.tanh(a::Measurement)
    aval = a.val
    return result(tanh(aval), abs2(sech(aval)), a)
end

# Other trig-related functions: sinpi, cospi, sinc, cosc

function Base.sinpi(a::Measurement)
    x = a.val
    return result(sinpi(x), cospi(x) * pi, a)
end

function Base.cospi(a::Measurement)
    x = a.val
    return result(cospi(x), -sinpi(x) * pi, a)
end

function Base.sinc(a::Measurement)
    x = a.val
    return result(sinc(x), cosc(x), a)
end

function Base.cosc(a::Measurement)
    x = a.val
    res = cosc(x)
    return result(res,
                  iszero(x) ? -oftype(x, pi) ^ 2 / 3 : -2 * res / x - sinc(x) * oftype(x, pi) ^ 2,
                  a)
end

# Inverse trig functions: acos, acosd, acosh, asin, asind, asinh, atan, atand, atanh,
#                         asec, acsc, acot, asech, acsch, acoth

function Base.acos(a::Measurement{T}) where {T<:AbstractFloat}
    aval = a.val
    return result(acos(aval), -inv(sqrt(one(T) - abs2(aval))), a)
end

function Base.acosd(a::Measurement{T}) where {T<:AbstractFloat}
    aval = a.val
    return result(acosd(aval), -rad2deg(inv(sqrt(one(T) - abs2(aval)))), a)
end

function Base.acosh(a::Measurement{T}) where {T<:AbstractFloat}
    aval = a.val
    return result(acosh(aval), inv(sqrt(abs2(aval) - one(T))), a)
end

function Base.asin(a::Measurement{T}) where {T<:AbstractFloat}
    aval = a.val
    return result(asin(aval), inv(sqrt(one(T) - abs2(aval))), a)
end

function Base.asind(a::Measurement{T}) where {T<:AbstractFloat}
    aval = a.val
    return result(asind(aval), rad2deg(inv(sqrt(one(T) - abs2(aval)))), a)
end

function Base.asinh(a::Measurement{T}) where {T<:AbstractFloat}
    aval = a.val
    return result(asinh(aval), inv(hypot(aval, one(T))), a)
end

function Base.atan(a::Measurement{T}) where {T<:AbstractFloat}
    aval = a.val
    return result(atan(aval), inv(abs2(aval) + one(T)), a)
end

function Base.atand(a::Measurement{T}) where {T<:AbstractFloat}
    aval = a.val
    return result(atand(aval), rad2deg(inv(abs2(aval) + one(T))), a)
end

function Base.atanh(a::Measurement{T}) where {T<:AbstractFloat}
    aval = a.val
    return result(atanh(aval), inv(one(T) - abs2(aval)), a)
end

function Base.atan(a::Measurement, b::Measurement)
    aval = a.val
    bval = b.val
    denom = abs2(aval) + abs2(bval)
    return result(atan(aval, bval),
                  (bval / denom, -aval / denom),
                  (a, b))
end

function Base.atan(a::Measurement, b::Real)
    x = a.val
    return result(atan(x, b), -b/(abs2(x) + abs2(b)), a)
end

function Base.atan(a::Real, b::Measurement)
    y = b.val
    return result(atan(a, y), -a/(abs2(a) + abs2(y)), b)
end

function Base.asec(a::Measurement)
    x = a.val
    return result(asec(x), 1 / (sqrt(x ^ 2 - 1) * abs(x)), a)
end

function Base.acsc(a::Measurement)
    x = a.val
    return result(acsc(x), -1 / (sqrt(x ^ 2 - 1) * abs(x)), a)
end

function Base.acot(a::Measurement)
    x = a.val
    return result(acot(x), -1 / (x ^ 2 + 1), a)
end

function Base.asech(a::Measurement)
    x = a.val
    return result(asech(x), -1 / (sqrt(1 - x ^ 2) * abs(x)), a)
end

function Base.acsch(a::Measurement)
    x = a.val
    return result(acsch(x), -1 / (sqrt(x ^ 2 + 1) * abs(x)), a)
end

function Base.acoth(a::Measurement)
    x = a.val
    return result(acoth(x), 1 / (1 - x ^ 2), a)
end

# Reciprocal trig functions: csc, cscd, csch, sec, secd, sech, cot, cotd, coth

function Base.csc(a::Measurement)
    aval = a.val
    val = csc(aval)
    return result(val, -val*cot(aval), a)
end

function Base.cscd(a::Measurement)
    aval = a.val
    val = cscd(aval)
    return result(val, -deg2rad(val*cotd(aval)), a)
end

function Base.csch(a::Measurement)
    aval = a.val
    val = csch(aval)
    return result(val, -val*coth(aval), a)
end

function Base.sec(a::Measurement)
    aval = a.val
    val = sec(aval)
    return result(val, val*tan(aval), a)
end

function Base.secd(a::Measurement)
    aval = a.val
    val = secd(aval)
    return result(val, deg2rad(val*tand(aval)), a)
end

function Base.sech(a::Measurement)
    aval = a.val
    val = sech(aval)
    return result(val, val*tanh(aval), a)
end

function Base.cot(a::Measurement)
    aval = a.val
    return result(cot(aval), -abs2(csc(aval)), a)
end

function Base.cotd(a::Measurement)
    aval = a.val
    return result(cotd(aval), -deg2rad(abs2(cscd(aval))), a)
end

function Base.coth(a::Measurement)
    aval = a.val
    return result(coth(aval), -abs2(csch(aval)), a)
end

### Exponential-related

# Exponentials: exp, expm1, exp10, frexp, ldexp

function Base.exp(a::Measurement)
    val = exp(a.val)
    return result(val, val, a)
end

function Base.expm1(a::Measurement)
    aval = a.val
    return result(expm1(aval), exp(aval), a)
end

function Base.exp10(a::Measurement{T}) where {T<:AbstractFloat}
    val = exp10(a.val)
    return result(val, log(T(10))*val, a)
end

function Base.frexp(a::Measurement)
    x, y = frexp(a.val)
    return (result(x, inv(exp2(y)), a), y)
end

Base.ldexp(a::Measurement{T}, e::Integer) where {T<:AbstractFloat} =
    result(ldexp(a.val, e), ldexp(one(T), e), a)

# Logarithms

function Base.log(a::Measurement, b::Measurement)
    aval = a.val
    bval = b.val
    val = log(aval, bval)
    loga = log(aval)
    return result(val, (-val / (aval * loga), 1 / (loga * bval)), (a, b))
end

function Base.log(a::Measurement) # Special case
    aval = a.val
    return result(log(aval), inv(aval), a)
end

function Base.log2(a::Measurement{T}) where {T<:AbstractFloat} # Special case
    x = a.val
    return result(log2(x), inv(log(T(2)) * x), a)
end

function Base.log10(a::Measurement{T}) where {T<:AbstractFloat} # Special case
    aval = a.val
    return result(log10(aval), inv(log(T(10)) * aval), a)
end

function Base.log1p(a::Measurement{T}) where {T<:AbstractFloat} # Special case
    aval = a.val
    return result(log1p(aval), inv(aval + one(T)), a)
end

Base.log(::Irrational{:ℯ}, a::Measurement) = log(a)

function Base.log(a::Real, b::Measurement)
    bval = b.val
    return result(log(a, bval), inv(log(a) * bval), b)
end

function Base.log(a::Measurement, b::Real)
    aval = a.val
    res = log(aval, b)
    return result(res, -res/(aval*log(aval)), a)
end

# Hypotenuse: hypot

function Base.hypot(a::Measurement, b::Measurement)
    aval = a.val
    bval = b.val
    val = hypot(aval, bval)
    return result(val,
                  (aval / val, bval / val),
                  (a, b))
end

function Base.hypot(a::Real, b::Measurement)
    bval = b.val
    res = hypot(a, bval)
    return result(res, bval / res, b)
end

function Base.hypot(a::Measurement, b::Real)
    aval = a.val
    res = hypot(aval, b)
    return result(res, aval / res, a)
end

# Square root: sqrt

function Base.sqrt(a::Measurement)
    val = sqrt(a.val)
    return result(val, inv(2 * val), a)
end

# Cube root: cbrt

function Base.cbrt(a::Measurement)
    aval = a.val
    val = cbrt(aval)
    return result(val, val / (3 * aval), a)
end

### Absolute value, sign and the likes

# Absolute value

function Base.abs(a::Measurement)
    aval = a.val
    return result(abs(aval), copysign(1, aval), a)
end

function Base.abs2(a::Measurement)
    x = a.val
    return result(abs2(x), 2*x, a)
end

# Sign: sign, copysign, flipsign

Base.sign(a::Measurement) = result(sign(a.val), 0, a)

Base.copysign(a::Measurement, b::Measurement) = ifelse(signbit(a)!=signbit(b), -a, a)
Base.copysign(a::Measurement, b::Real)        = ifelse(signbit(a)!=signbit(b), -a, a)
Base.flipsign(a::Measurement, b::Measurement) = ifelse(signbit(b), -a, a)
Base.flipsign(a::Measurement, b::Real)        = ifelse(signbit(b), -a, a)
for T in (Signed, Rational, Float32, Float64, Real)
    @eval Base.copysign(a::$T, b::Measurement) = copysign(a, b.val)
    @eval Base.flipsign(a::$T, b::Measurement) = flipsign(a, b.val)
end

### Modulo

# Use definition of "mod" function:
# http://docs.julialang.org/en/stable/manual/mathematical-operations/#division-functions
Base.mod(a::Measurement, b::Measurement) = a - fld(a, b)*b
Base.mod(a::Measurement, b::Real) = result(mod(a.val, b), 1, a)
Base.mod(a::Real, b::Measurement) = result(mod(a, b.val), -fld(a, b.val), b)

# Use definition of "rem" function:
# https://docs.julialang.org/en/latest/manual/mathematical-operations/#Division-functions-1
Base.rem(a::Measurement, b::Measurement) = a - div(a, b)*b
Base.rem(a::Measurement, b::Real) = result(rem(a.val, b), 1, a)
Base.rem(a::Real, b::Measurement) = result(rem(a, b.val), -div(a, b.val), b)
Base.rem(a::Measurement, b::Union{Measurement,Float64}, ::RoundingMode{:Nearest}) =
    a - b * round(a / b, RoundNearest)
Base.rem(a::Float64, b::Measurement, ::RoundingMode{:Nearest}) =
    a - b * round(a / b, RoundNearest)

Base.mod2pi(a::Measurement) = result(mod2pi(a.val), 1, a)
Base.rem2pi(a::Measurement, r::RoundingMode) = result(rem2pi(a.val, r), 1, a)

### Machine precision

Base.eps(::Type{Measurement{T}}) where {T<:AbstractFloat} = eps(T)
Base.eps(a::Measurement) = eps(a.val)

Base.nextfloat(a::Measurement) = result(nextfloat(a.val), 1, a)
Base.nextfloat(a::Measurement, n::Integer) = result(nextfloat(a.val, n), 1, a)

Base.maxintfloat(::Type{Measurement{T}}) where {T<:AbstractFloat} = maxintfloat(T)

Base.floatmin(::Type{Measurement{T}}) where {T<:AbstractFloat} = floatmin(T) ± zero(T)
Base.floatmax(::Type{Measurement{T}}) where {T<:AbstractFloat} = floatmax(T) ± zero(T)

Base.typemax(::Type{Measurement{T}}) where {T<:AbstractFloat} = typemax(T)

### Rounding

Base.round(a::Measurement, r::RoundingMode=RoundNearest; kwargs...) =
    measurement(round(value(a), r; kwargs...), round(uncertainty(a); kwargs...))
Base.round(::Type{T}, a::Measurement, r::RoundingMode=RoundNearest) where {T<:Integer} =
    round(T, a.val, r)
Base.floor(a::Measurement) = measurement(floor(a.val))
Base.floor(::Type{T}, a::Measurement) where {T<:Integer} = floor(T, a.val)
Base.ceil(a::Measurement) = measurement(ceil(a.val))
Base.ceil(::Type{T}, a::Measurement) where {T<:Integer} = ceil(Integer, a.val)
Base.trunc(a::Measurement) = measurement(trunc(a.val))
Base.trunc(::Type{T}, a::Measurement) where {T<:Integer} = trunc(T, a.val)

# Widening

Base.widen(::Type{Measurement{T}}) where {T<:AbstractFloat} = Measurement{widen(T)}

# To big float

Base.big(::Type{Measurement}) = Measurement{BigFloat}
Base.big(::Type{Measurement{T}}) where {T<:AbstractFloat} = Measurement{BigFloat}
Base.big(x::Measurement{<:AbstractFloat}) = convert(Measurement{BigFloat}, x)
Base.big(x::Complex{<:Measurement}) = convert(Complex{Measurement{BigFloat}}, x)

# Sum and prod

# This definition is not strictly needed, because `sum' works out-of-the-box
# with Measurement type, but this makes the function linear instead of quadratic
# in the number of arguments, but `result' is quadratic in the number of
# arguments, so in the end the function goes from cubic to quadratic.  Still not
# ideal, but this is an improvement.
Base.sum(a::AbstractArray{<:Measurement}) = result(sum(value.(a)), ones(length(a)), a)

# Same as above.  I'm not particularly proud of how the derivatives are
# computed, but something like this is needed in order to avoid errors with null
# nominal values: you may think to x ./ prod(x), but that would fail if one or
# more elements are zero.
function Base.prod(a::AbstractArray{<:Measurement})
    x = value.(a)
    return result(prod(x),
                  [prod(deleteat!(copy(x), i)) for i in eachindex(x)],
                  a)
end


text:
module TestMatrices
using Test
using SparseArrays

using Catlab.CategoricalAlgebra.Matrices

# Dense matrices
################

const IntMat = Matrix{Int}

A, B, C = [1 2 3; 4 5 6], [1 -1; -1 2], [0 1; 1 0]
@test dom(A) == MatrixDom{IntMat}(3)
@test codom(A) == MatrixDom{IntMat}(2)
@test A⋅B == B*A
@test id(dom(A))⋅A == A
@test A⋅id(codom(A)) == A

@test A⊕(B⊕C) == (A⊕B)⊕C
@test swap(dom(A),dom(B))⋅(B⊕A)⋅swap(codom(B),codom(A)) == A⊕B
@test pair(B,C) == mcopy(dom(B))⋅(B⊕C)
@test copair(A,B) == (A⊕B)⋅plus(codom(A))

@test A⊗(B⊗C) == (A⊗B)⊗C
@test braid(dom(A),dom(B))⋅(B⊗A)⋅braid(codom(B),codom(A)) == A⊗B

# Sparse matrices
#################

const SparseIntMat = SparseMatrixCSC{Int,Int}

A, B, C = sparse(A), sparse(B), sparse(C)
@test dom(A) == MatrixDom{SparseIntMat}(3)
@test codom(A) == MatrixDom{SparseIntMat}(2)
@test A⋅B isa SparseIntMat
@test A⊕B isa SparseIntMat
@test A⊗B isa SparseIntMat
@test A⋅B == B*A
@test id(dom(A))⋅A == A
@test A⋅id(codom(A)) == A

end


text:
using Random
using Images, CoordinateTransformations
using ImageView

dic = Dict()
arr = ones(3,2,2)/2.
@show arr
img = colorview(RGB, arr)
dic[1]=img
@show img
chview = channelview(img)
chview[1]=0
@show arr
@show dic[1]

max_stars_repo_name: terasakisatoshi/juliaExer

max_issues_repo_name: terasakisatoshi/juliaExer

max_forks_repo_name: terasakisatoshi/juliaExer


text:
# Determine the local system's time zone
# Based upon Python's tzlocal https://pypi.python.org/pypi/tzlocal
import Compat: @static, is_apple, is_unix, is_windows, readstring
using Mocking

"""
    localzone() -> TimeZone

Returns a `TimeZone` object that is equivalent to the system's current time zone.
"""
function localzone()
    @static if is_apple()
        name = @mock readstring(`systemsetup -gettimezone`)  # Appears to only work as root
        if contains(name, "Time Zone: ")
            name = strip(replace(name, "Time Zone: ", ""))
        else
            # link will be something like /usr/share/zoneinfo/Europe/Warsaw
            name = @mock readlink("/etc/localtime")
            name = match(r"(?<=zoneinfo/).*$", name).match
        end
        return TimeZone(name)
    elseif is_unix()
        name = ""
        validnames = timezone_names()

        # Try getting the time zone from the "TZ" environment variable
        # http://linux.die.net/man/3/tzset
        if haskey(ENV, "TZ")
            name = ENV["TZ"]
            startswith(name, ':') || error("Currently only support filespec for TZ variable")
            name = name[2:end]

            if startswith(name, '/')
                return @mock open(name) do f
                    read_tzfile(f, "local")
                end
            else
                # Relative name matches pre-compiled time zone name
                name in validnames && return TimeZone(name)

                # The system time zone directory used depends on the (g)libc version
                tzdirs = ["/usr/lib/zoneinfo", "/usr/share/zoneinfo"]
                haskey(ENV, "TZDIR") && unshift!(tzdirs, ENV["TZDIR"])

                for dir in tzdirs
                    filepath = joinpath(dir, name)
                    (@mock isfile(filepath)) || continue
                    return @mock open(filepath) do f
                        read_tzfile(f, name)
                    end
                end

                throw(SystemError("unable to locate tzfile: $name"))
            end
        end

        # Look for distribution specific configuration files that contain the time zone name.

        filename = "/etc/timezone"
        if @mock isfile(filename)
            @mock open(filename) do file
                name = readstring(file)

                # Get rid of host definitions and comments:
                name = strip(replace(name, r"#.*", ""))
                name = replace(name, ' ', '_')
            end

            name in validnames && return TimeZone(name)
        end

        # CentOS has a ZONE setting in /etc/sysconfig/clock,
        # OpenSUSE has a TIMEZONE setting in /etc/sysconfig/clock and
        # Gentoo has a TIMEZONE setting in /etc/conf.d/clock

        zone_re = r"(TIME)?ZONE\s*=\s*\"(?<name>.*?)\""
        for filepath in ("/etc/sysconfig/clock", "/etc/conf.d/clock")
            (@mock isfile(filepath)) || continue
            @mock open(filepath) do file
                for line in readlines(file)
                    matched = match(zone_re, line)
                    if matched != nothing
                        name = matched["name"]
                        name = replace(name, ' ', '_')
                        break
                    end
                end
            end

            name in validnames && return TimeZone(name)
        end

        # systemd distributions use symlinks that include the zone name,
        # see manpage of localtime(5) and timedatectl(1)
        link = "/etc/localtime"
        if @mock islink(link)
            filepath = @mock readlink(link)
            start = search(filepath, '/')

            while start != 0
                name = filepath[(start + 1):end]
                name in validnames && return TimeZone(name)
                start = search(filepath, '/', start + 1)
            end
        end

        # No explicit setting existed. Use localtime
        for filepath in ("/etc/localtime", "/usr/local/etc/localtime")
            (@mock isfile(filepath)) || continue
            return @mock open(filepath) do f
                read_tzfile(f, "local")
            end
        end
    elseif is_windows()
        isfile(WIN_TRANSLATION_FILE) || error("Missing Windows to POSIX time zone translation ",
            "file. Try running Pkg.build(\"TimeZones\")")

        translation = open(WIN_TRANSLATION_FILE, "r") do fp
            deserialize(fp)
        end

        # Windows powershell should be available on Windows 7 and above
        win_name = strip(@mock readstring(`powershell -Command "[TimeZoneInfo]::Local.Id"`))
        if haskey(translation, win_name)
            posix_name = translation[win_name]

            # Translation dict includes Etc time zones which we currently are not supporting
            # since they are deemed historical. To ensure compatibility with the translation
            # dict we will manually convert these fixed time zones.
            if startswith(posix_name, "Etc/GMT")
                name = replace(posix_name, r"Etc/GMT0?", "UTC")

                # Note: Etc/GMT[+-] are reversed compared to UTC[+-]
                if contains(name, "+")
                    name = replace(name, "+", "-")
                else
                    name = replace(name, "-", "+")
                end

                return FixedTimeZone(name)
            else
                return TimeZone(posix_name)
            end
        else
            error("unable to translate to POSIX time zone name from: \"$win_name\"")
        end
    end

    error("Failed to find local time zone")
end


text:
export iauSp00
"""
The TIO locator s', positioning the Terrestrial Intermediate Origin
on the equator of the Celestial Intermediate Pole.

This function is part of the International Astronomical Union's
SOFA (Standards Of Fundamental Astronomy) software collection.

Status:  canonical model.

Given:
   date1,date2  double    TT as a 2-part Julian Date (Note 1)

Returned (function value):
               double    the TIO locator s' in radians (Note 2)

Notes:

   1. The TT date date1+date2 is a Julian Date, apportioned in any
   convenient way between the two arguments.  For example,
   JD(TT)=2450123.7 could be expressed in any of these ways,
   among others:

         date1          date2

      2450123.7           0.0       (JD method)
      2451545.0       -1421.3       (J2000 method)
      2400000.5       50123.2       (MJD method)
      2450123.5           0.2       (date & time method)

   The JD method is the most natural and convenient to use in
   cases where the loss of several decimal digits of resolution
   is acceptable.  The J2000 method is best matched to the way
   the argument is handled internally and will deliver the
   optimum resolution.  The MJD method and the date & time methods
   are both good compromises between resolution and convenience.

   2. The TIO locator s' is obtained from polar motion observations by
   numerical integration, and so is in essence unpredictable.
   However, it is dominated by a secular drift of about
   47 microarcseconds per century, which is the approximation
   evaluated by the present function.

Reference:

   McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003),
   IERS Technical Note No. 32, BKG (2004)

This revision:  2013 June 18

SOFA release 2018-01-30

Copyright (C) 2018 IAU SOFA Board.  See notes at end.
"""

# double iauSp00(double date1, double date2)

function iauSp00(date1::Real, date2::Real)
   return ccall((:iauSp00, libsofa_c), Cdouble,
               (Cdouble, Cdouble),
               convert(Float64, date1), convert(Float64, date2))
end


text:
using Test, FinancialMonteCarlo,Statistics;
@show "Black Scholes Model"
S0=100.0;
K=100.0;
r=[0.00,0.02];
T=1.0;
d=0.01;
D=90.0;

Nsim=100000;
Nstep=30;
sigma=0.2;
mc=MonteCarloConfiguration(Nsim,Nstep);
toll=0.8

rfCurve=ZeroRate(r,T);

FwdData=Forward(T)
EUData=EuropeanOption(T,K)

S=S0.*exp.(simulate(BrownianMotion(sigma,-sigma*sigma*0.5,Underlying(0.0,0.0)),rfCurve,mc,T));

Payoff=exp(-integral(rfCurve,T))*S[:,end];
Price=mean(Payoff)


text:
### Somatório de idade
# Uma criança quer saber qual é a soma de todas as idades que ela já teve.
# Elaborar programa que lê  uma idade qualquer e
# responde rapidamente a essa pergunta
# ( fórmula para calcular a soma dos N primeiros números inteiros: N (N+1)/2 ).

print("inform a age ")
age = parse(UInt8, readline())
age = age * ((age+1)/2)

print(UInt8(age),"\n")

max_stars_repo_name: codinginbrazil/functional-programming

max_issues_repo_name: codinginbrazil/functional-programming

max_forks_repo_name: codinginbrazil/functional-programming


text:
"""
    FlipSignBridge{S1, S2, F}

Bridge a `G`-in-`S1` constraint into an `F`-in-`S2` constraint by multiplying
the function by `-1` and taking the point reflection of the set across the
origin. The flipped `F`-in-`S` constraint is stored in the `flipped_constraint`
field by convention.
"""
abstract type FlipSignBridge{
    S1<:MOI.AbstractSet, S2<:MOI.AbstractSet,
    F<:MOI.AbstractFunction} <: AbstractBridge end

function MOI.supports_constraint(::Type{<:FlipSignBridge{S1}},
                                 ::Type{<:MOI.AbstractScalarFunction},
                                 ::Type{S1}) where {S1<:MOI.AbstractSet}
    return true
end
function added_constraint_types(
    ::Type{<:FlipSignBridge{S1, S2, F}}) where {S1, S2, F}
    return [(F, S2)]
end

# Attributes, Bridge acting as an model
function MOI.get(::FlipSignBridge{S1, S2, F},
                 ::MOI.NumberOfConstraints{F, S2}) where {S1, S2, F}
    return 1
end
function MOI.get(bridge::FlipSignBridge{S1, S2, F},
                 ::MOI.ListOfConstraintIndices{F, S2}) where {S1, S2, F}
    return [bridge.flipped_constraint]
end

# References
function MOI.delete(model::MOI.ModelLike, bridge::FlipSignBridge)
    MOI.delete(model, bridge.flipped_constraint)
end

# Attributes, Bridge acting as a constraint
function MOI.get(model::MOI.ModelLike,
                 attr::Union{MOI.ConstraintPrimal, MOI.ConstraintDual},
                 bridge::FlipSignBridge)
    return -MOI.get(model, attr, bridge.flipped_constraint)
end

function MOI.modify(model::MOI.ModelLike, bridge::FlipSignBridge,
                    change::MOI.ScalarCoefficientChange)
    MOI.modify(
        model, bridge.flipped_constraint,
        MOI.ScalarCoefficientChange(change.variable, -change.new_coefficient))
end

function MOI.modify(model::MOI.ModelLike, bridge::FlipSignBridge,
                    change::MOI.MultirowChange{T}) where T
    new_coefficients = Tuple{Int64, T}[
        (index, -coef) for (index, coef) in change.new_coefficients]
    MOI.modify(model, bridge.flipped_constraint,
               MOI.MultirowChange(change.variable,
                                  new_coefficients))
end
function MOI.modify(model::MOI.ModelLike, bridge::FlipSignBridge,
                    change::MOI.VectorConstantChange)
    MOI.modify(model, bridge.flipped_constraint,
               MOI.VectorConstantChange(-change.new_constant))
end

"""
    GreaterToLessBridge{T, F<:MOI.AbstractScalarFunction} <:
        FlipSignBridge{MOI.GreaterThan{T}, MOI.LessThan{T}, F}

Transforms a `G`-in-`GreaterThan{T}` constraint into an `F`-in-`LessThan{T}`
constraint.
"""
struct GreaterToLessBridge{T, F<:MOI.AbstractScalarFunction} <:
    FlipSignBridge{MOI.GreaterThan{T}, MOI.LessThan{T}, F}
    flipped_constraint::CI{F, MOI.LessThan{T}}
end
function GreaterToLessBridge{T, F}(model::MOI.ModelLike,
                                   g::MOI.AbstractScalarFunction,
                                   s::MOI.GreaterThan) where {T, F}
    f = MOIU.operate(-, T, g)
    flipped_constraint = MOI.add_constraint(model, f, MOI.LessThan(-s.lower))
    return GreaterToLessBridge{T, F}(flipped_constraint)
end
function concrete_bridge_type(::Type{<:GreaterToLessBridge{T}},
                              G::Type{<:MOI.AbstractScalarFunction},
                              ::Type{MOI.GreaterThan{T}}) where T
    F = MOIU.promote_operation(-, T, G)
    return GreaterToLessBridge{T, F}
end
function MOI.set(model::MOI.ModelLike, attr::MOI.ConstraintSet,
                 bridge::GreaterToLessBridge, new_set::MOI.GreaterThan)
    MOI.set(model, attr, bridge.flipped_constraint,
            MOI.LessThan(-new_set.lower))
end

"""
    LessToGreaterBridge{T, F<:MOI.AbstractScalarFunction} <:
        FlipSignBridge{MOI.LessThan{T}, MOI.GreaterThan{T}, F}

Transforms a `G`-in-`LessThan{T}` constraint into an `F`-in-`GreaterThan{T}`
constraint.
"""
struct LessToGreaterBridge{T, F<:MOI.AbstractScalarFunction} <:
    FlipSignBridge{MOI.LessThan{T}, MOI.GreaterThan{T}, F}
    flipped_constraint::CI{F, MOI.GreaterThan{T}}
end
function LessToGreaterBridge{T, F}(model::MOI.ModelLike,
                                   g::MOI.AbstractScalarFunction,
                                   s::MOI.LessThan) where {T, F}
    f = MOIU.operate(-, T, g)
    flipped_constraint = MOI.add_constraint(model, f, MOI.GreaterThan(-s.upper))
    return LessToGreaterBridge{T, F}(flipped_constraint)
end
function concrete_bridge_type(::Type{<:LessToGreaterBridge{T}},
                              G::Type{<:MOI.AbstractScalarFunction},
                              ::Type{MOI.LessThan{T}}) where T
    F = MOIU.promote_operation(-, T, G)
    return LessToGreaterBridge{T, F}
end
function MOI.set(model::MOI.ModelLike, attr::MOI.ConstraintSet,
                 bridge::LessToGreaterBridge, new_set::MOI.LessThan)
    MOI.set(model, attr, bridge.flipped_constraint,
            MOI.GreaterThan(-new_set.upper))
end

"""
    NonnegToNonposBridge{T, F<:MOI.AbstractVectorFunction}

Transforms a `G`-in-`Nonnegatives` constraint into a `F`-in-`Nonpositives`
constraint.
"""
struct NonnegToNonposBridge{T, F<:MOI.AbstractVectorFunction} <:
    FlipSignBridge{MOI.Nonnegatives, MOI.Nonpositives, F}
    flipped_constraint::CI{F, MOI.Nonpositives}
end
function NonnegToNonposBridge{T, F}(model::MOI.ModelLike,
                                    g::MOI.AbstractVectorFunction,
                                    s::MOI.Nonnegatives) where {T, F}
    f = MOIU.operate(-, T, g)
    flipped_constraint = MOI.add_constraint(model, f,
                                            MOI.Nonpositives(s.dimension))
    return NonnegToNonposBridge{T, F}(flipped_constraint)
end
function concrete_bridge_type(::Type{<:NonnegToNonposBridge{T}},
                              G::Type{<:MOI.AbstractVectorFunction},
                              ::Type{MOI.Nonnegatives}) where T
    F = MOIU.promote_operation(-, T, G)
    return NonnegToNonposBridge{T, F}
end

"""
    NonposToNonnegBridge{T, F<:MOI.AbstractVectorFunction}

Transforms a `G`-in-`Nonpositives` constraint into a `F`-in-`Nonnegatives`
constraint.
"""
struct NonposToNonnegBridge{T, F<:MOI.AbstractVectorFunction} <:
    FlipSignBridge{MOI.Nonpositives, MOI.Nonnegatives, F}
    flipped_constraint::CI{F, MOI.Nonnegatives}
end
function NonposToNonnegBridge{T, F}(model::MOI.ModelLike,
                                    g::MOI.AbstractVectorFunction,
                                    s::MOI.Nonpositives) where {T, F}
    f = MOIU.operate(-, T, g)
    flipped_constraint = MOI.add_constraint(model, f,
                                            MOI.Nonnegatives(s.dimension))
    return NonposToNonnegBridge{T, F}(flipped_constraint)
end
function concrete_bridge_type(::Type{<:NonposToNonnegBridge{T}},
                              G::Type{<:MOI.AbstractVectorFunction},
                              ::Type{MOI.Nonpositives}) where T
    F = MOIU.promote_operation(-, T, G)
    return NonposToNonnegBridge{T, F}
end


text:
#arc diagrams

# plots
#import Pkg
#Pkg.add("LinearAlgebra")
#Pkg.add("GraphRecipes")
#Pkg.add("SparseArrays")

using LinearAlgebra
using SparseArrays
using GraphRecipes
using Plots

adjmat = Symmetric(sparse(rand(0:1,8,8)))

plot(
    graphplot(adjmat,
              method=:arcdiagram,
              markersize=0.5,
              linecolor=:black,
              markercolor=:black)
    )

max_stars_repo_name: Ellon-M/visualizations

max_issues_repo_name: Ellon-M/visualizations

max_forks_repo_name: Ellon-M/visualizations


text:
using DifferentialEquations
using DiffEqSensitivity
using Random
using Distributions
using Turing
using DataFrames
using StatsPlots

function sir_ode!(du,u,p,t)
    (S,I,R,C) = u
    (β,c,γ) = p
    N = S+I+R
    infection = β*c*I/N*S
    recovery = γ*I
    @inbounds begin
        du[1] = -infection
        du[2] = infection - recovery
        du[3] = recovery
        du[4] = infection
    end
    nothing
end;

tmax = 40.0
tspan = (0.0,tmax)
obstimes = 1.0:1.0:tmax
u0 = [990.0,10.0,0.0,0.0] # S,I.R,C
p = [0.05,10.0,0.25]; # β,c,γ

prob_ode = ODEProblem(sir_ode!,u0,tspan,p);

sol_ode = solve(prob_ode,
            Tsit5(),
            saveat = 1.0);

C = Array(sol_ode)[4,:] # Cumulative cases
X = C[2:end] - C[1:(end-1)];

Random.seed!(1234)
Y = rand.(Poisson.(X));

bar(obstimes,Y,legend=false)
plot!(obstimes,X,legend=false)

@model bayes_sir(y) = begin
  # Calculate number of timepoints
  l = length(y)
  i₀  ~ Uniform(0.0,1.0)
  β ~ Uniform(0.0,1.0)
  I = i₀*1000.0
  u0=[1000.0-I,I,0.0,0.0]
  p=[β,10.0,0.25]
  tspan = (0.0,float(l))
  prob = ODEProblem(sir_ode!,
          u0,
          tspan,
          p)
  sol = solve(prob,
              Tsit5(),
              saveat = 1.0)
  sol_C = Array(sol)[4,:] # Cumulative cases
  sol_X = sol_C[2:end] - sol_C[1:(end-1)]
  l = length(y)
  for i in 1:l
    y[i] ~ Poisson(sol_X[i])
  end
end;

ode_nuts = sample(bayes_sir(Y),NUTS(0.65),10000);

describe(ode_nuts)

plot(ode_nuts)

posterior = DataFrame(ode_nuts);

histogram2d(posterior[!,:β],posterior[!,:i₀],
                bins=80,
                xlabel="β",
                ylab="i₀",
                ylim=[0.006,0.016],
                xlim=[0.045,0.055],
                legend=false)
plot!([0.05,0.05],[0.0,0.01])
plot!([0.0,0.05],[0.01,0.01])

function predict(y,chain)
    # Length of data
    l = length(y)
    # Length of chain
    m = length(chain)
    # Choose random
    idx = sample(1:m)
    i₀ = chain[:i₀][idx]
    β = chain[:β][idx]
    I = i₀*1000.0
    u0=[1000.0-I,I,0.0,0.0]
    p=[β,10.0,0.25]
    tspan = (0.0,float(l))
    prob = ODEProblem(sir_ode!,
            u0,
            tspan,
            p)
    sol = solve(prob,
                Tsit5(),
                saveat = 1.0)
    out = Array(sol)
    sol_X = [0.0; out[4,2:end] - out[4,1:(end-1)]]
    hcat(sol_ode.t,out',sol_X)
end;

Xp = []
for i in 1:10
    pred = predict(Y,ode_nuts)
    push!(Xp,pred[2:end,6])
end

scatter(obstimes,Y,legend=false)
plot!(obstimes,Xp,legend=false)


text:
using ImageFiltering, ImageCore, OffsetArrays, Logging, ImageMetadata, Test
import AxisArrays
using AxisArrays: AxisArray, Axis

@testset "basic" begin
    v = OffsetArray([1,2,3], -1:1)
    @test reflect(v) == OffsetArray([3,2,1], -1:1)

    rv = KernelFactors.ReshapedOneD{3,1}(v)
    @test ndims(rv) == 3
    @test rv[0] == 2
    @test rv == KernelFactors.ReshapedOneD{3,1}(Float64.(v))
    @test rv != KernelFactors.ReshapedOneD{3,1}([1,2,3])
    @test reshape(v, 0:0, -1:1, 0:0) == rv
    k3 = (KernelFactors.ReshapedOneD{3,0}(v), rv, KernelFactors.ReshapedOneD{3,2}(v))
    @test ImageFiltering.isseparable(k3)
    @test .*(k3...) == v .* reshape(v, 0:0, -1:1) .* reshape(v, 0:0, 0:0, -1:1)

    kern = OffsetArray([1 4 7; 2 5 8; 3 6 9], 0:2, -1:1)
    @test reflect(kern) == OffsetArray([9 6 3; 8 5 2; 7 4 1], -2:0, -1:1)

    @test eltype(KernelFactors.IIRGaussian(3)) == Float64
    @test eltype(KernelFactors.IIRGaussian(Float32, 3)) == Float32
    kern = KernelFactors.IIRGaussian([1,2.0f0])
    @test isa(kern, Tuple{KernelFactors.ReshapedOneD{Float32},KernelFactors.ReshapedOneD{Float32}})
    @test isa(ImageFiltering.filter_algorithm([1],[1],kern), Algorithm.IIR)
    @test isa(ImageFiltering.filter_algorithm([1],[1],(kern...,Kernel.Laplacian())), Algorithm.Mixed)
    kern = KernelFactors.IIRGaussian(Float64, [1,2.0f0])
    @test isa(kern, Tuple{KernelFactors.ReshapedOneD{Float64},KernelFactors.ReshapedOneD{Float64}})

    @test ndims(Pad(:replicate, (3,3))) == 2

    @test KernelFactors.kernelfactors(([0,3], [1,7]))  == (reshape([0,3], 1:2, 0:0), reshape([1,7], 0:0, 1:2))
    @test KernelFactors.kernelfactors(([0,3], [1,7]')) == (reshape([0,3], 2, 1), reshape([1,7], 1, 2))

    tiles = ImageFiltering.tile_allocate(Float32, (rand(3),rand(3)'))
    @test isa(tiles, Vector{Matrix{Float32}})
    @test length(tiles) == Threads.nthreads()

    @test length(ImageFiltering.safetail(CartesianIndices(()))) == 1
    @test ImageFiltering.safetail(CartesianIndex(())) == CartesianIndex(())
    @test length(ImageFiltering.safehead(CartesianIndices(()))) == 1
    @test ImageFiltering.safehead(CartesianIndex(())) == CartesianIndex(())

    # Warnings
    fname = tempname()
    open(fname, "w") do f
        logger = SimpleLogger(f)
        with_logger(logger) do
            KernelFactors.IIRGaussian(0.5, emit_warning=false)
        end
    end
    @test isempty(chomp(read(fname, String)))
    open(fname, "w") do f
        logger = SimpleLogger(f)
        with_logger(logger) do
            KernelFactors.IIRGaussian(0.5)
        end
    end
    @test occursin("too small for accuracy", read(fname, String))
    rm(fname)
end

@testset "centered" begin
    check_range(r, f, l) = (@test first(r) == f; @test last(r) == l)
    check_range_axes(r, f, l) = check_range(axes(r)[1], f, l)

    check_range(axes(centered(1:3))[1], -1, 1)
    a = AxisArray(rand(3, 3), Axis{:y}(0.1:0.1:0.3), Axis{:x}(1:3))
    ca = centered(a)
    axs = axes(ca)
    check_range(axs[1], -1, 1)
    check_range(axs[2], -1, 1)
    axs = AxisArrays.axes(ca)
    check_range(axs[1].val, 0.1, 0.3)
    check_range(axs[2].val, 1, 3)
    check_range_axes(axs[1].val, -1, 1)
    check_range_axes(axs[1].val, -1, 1)
    am = ImageMeta(a; prop1="simple")
    ca = centered(am)
    axs = axes(ca)
    check_range(axs[1], -1, 1)
    check_range(axs[2], -1, 1)
    axs = AxisArrays.axes(ca)
    check_range(axs[1].val, 0.1, 0.3)
    check_range(axs[2].val, 1, 3)
    check_range_axes(axs[1].val, -1, 1)
    check_range_axes(axs[1].val, -1, 1)

    a = rand(3, 3)
    ca = centered(a)
    cca = centered(ca)
    c_view_a = centered(view(ca, :, 0:0))
    c_slice_a = centered(a[:, 2:2])
    @test a[2, 2] == ca[0, 0]
    @test a[2, 2] == cca[0, 0]
    @test a[2, 2] == c_view_a[0, 0]
    @test a[2, 2] == c_slice_a[0, 0]
end

@testset "freqkernel/spacekernel" begin
    k = centered(reshape([1,-1], 2, 1))
    kfft = freqkernel(k, (31, 31))
    @test size(kfft) == (31, 31)
    k2 = real.(spacekernel(kfft, axes(k)))
    @test k2 ≈ k
    kfft = freqkernel(k, (31, 31); rfft=true)
    @test size(kfft) == (16, 31)
    k2 = spacekernel(kfft, axes(k); rfftsz=31)
    @test k2 ≈ k

    for T in (Float64, Float32, Float16, N0f16)
        k = T.(Kernel.gaussian(3))
        kfft = freqkernel(k)
        @test size(kfft) == size(k)
        k2 = real.(spacekernel(kfft))
        @test T.(k2) ≈ k

        kfft = freqkernel(k; rfft=true)
        k2 = spacekernel(kfft; rfftsz=size(k, 1))
        @test T.(k2) ≈ k
    end
end


text:
using Distributed
addprocs(16)
@everywhere include("plot_Jray.jl")
compute_rayleigh_criterion(1120)

max_stars_repo_name: nishaChandramoorthy/energies

max_issues_repo_name: nishaChandramoorthy/energies

max_forks_repo_name: nishaChandramoorthy/energies


text:
using Statistics, FixedEffectModels, DataFrames

struct EqmData
    Y::Array{Float64}
    D::Array{Float64}
    S::Array{Float64}
    Z::Array{Float64}
    W::Array{Bool, 1}
    U::Array{Float64}
    p::Float64
end

function EqmData(Y, D, S, W, U; p=0.5)
    return EqmData(Y, D, S, D .- S, W, U, p)
end

struct Partials
    yp1::Float64
    yp0::Float64
    dp1::Float64
    dp0::Float64
    sp1::Float64
    sp0::Float64
    μzp::Float64
end

function Partials(Y, D, S, W, U; p=0.5)
    ∂yp1, ∂yp0 = deriv_estimate(Y, W, U)
    ∂dp1, ∂dp0 = deriv_estimate(D, W, U)
    ∂sp1, ∂sp0 = deriv_estimate(S, W, U)
    ∂μzp = p*(∂dp1 - ∂sp1) + (1 - p)*(∂dp0 - ∂sp0)
    return Partials(∂yp1, ∂yp0, ∂dp1, ∂dp0, ∂sp1, ∂sp0, ∂μzp)
end

struct DiffMeans
    y::Float64
    d::Float64
    s::Float64
    z::Float64
end

function DiffMeans(Y, D, S, W; p=0.5)
    y = horvitz_thompson(Y, W)
    d = horvitz_thompson(D, W)
    s = horvitz_thompson(S, W)
    DiffMeans(y, d, s, d - s)
end

function horvitz_thompson(Y, W; p=0.5)
    data = (W./p - (1 .-W)./(1 - p)).*Y
    return mean(data)
end

function deriv_estimate(Y, W, U)
    ∂y1 = cov(Y[W.==1], U[W.==1])/var(U[W.==1])
    ∂y0 = cov(Y[W.==0], U[W.==0])/var(U[W.==0])
    return (∂y1, ∂y0)
end

function varianceDirect(Y, Z, W, ∂::Partials; p=0.5)
    a = (∂.yp1 - ∂.yp0)
    H = ∂.μzp
    A = -1*p .* W.*a./H.*Z + ((1 - p)*a/H) .* (1 .- W).*Z
    data = (W./p - (1 .-W)./(1 - p)).*(Y.+A)
    return std(data)/sqrt(length(Y))
end

function estimate_τ(data::EqmData)
    ∂ = Partials(data.Y, data.D, data.S, data.W, data.U)
    Δ = DiffMeans(data.Y, data.D, data.S, data.W)
    τD = Δ.y
    #println((data.p*∂.yp1 + (1- data.p)*∂.yp0))
    τI = -1.0*(data.p*∂.yp1 + (1- data.p)*∂.yp0)*Δ.z/∂.μzp
    σD = varianceDirect(data.Y, data.Z, data.W, ∂)
    σI = varianceIndirect(data.Y, data.Z, data.W, data.U, ∂, Δ)
    return (τD = τD, σD = σD, τI = τI, σI = σI)
end

function ivCOV(Y, X, Z)
    df = DataFrame(Y=Y, X=X, Z=Z)
    regmodel = reg(df, @formula(Y ~ (X ~ Z)), Vcov.robust())
    return (coef(regmodel)[2], stderror(regmodel)[2], coef(regmodel)[1])
end
function varianceIndirect(Y, Z, W, U, ∂::Partials, Δ::DiffMeans; p = 0.5)
    dpstar = -1.0*Δ.z
    coefIV, covIV, int = ivCOV(Y, Z, U)
    a = (∂.yp1 - ∂.yp0)
    H = ∂.μzp
    A = -1*p .* W.*a./H.*Z + ((1 - p)*a/H) .* (1 .- W).*Z
    data = (W./p - (1 .-W)./(1 - p)).*(Z.+A)
    covZ = var(data)/length(data)
    #σ2 = dpstar*(covIV^2*dpstar) + coefIV^2*covZ
    #println("Cov direct:", sqrt(dpstar*(covIV^2*dpstar)))
    ehat = Y .- Z.* coefIV .- int
    estimatehat = mean(ehat.^2)*(Δ.z/∂.μzp)^2/(length(Y)*U[1]^2)
    σ2 = estimatehat + coefIV^2*covZ
    return sqrt(σ2)#/abs(U[1])
end

max_stars_repo_name: evanmunro/market-equilibrium

max_issues_repo_name: evanmunro/market-equilibrium

max_forks_repo_name: evanmunro/market-equilibrium