text:
function RTs_SS(a,b,distri,cohmeasure)
    # a quantity to mean
    # b quantity to test
    # Return the mean of RTs and STD
    # Sucess after Sucess
    output = Float64[]
    for i=1:length(a)-1
        x = b[i]
        y=b[i+1]
        x > 0 && y >0 && distri[i+1]*distri[i] > 0 && distri[i+1] == cohmeasure || continue
        push!(output,a[i+1])
    end
    return mean(output),std(output),length(output)
end

function performance(results)
    # compute the performance of the model
    # need the list of of results for a specific coherence
    return sum(results)/length(results)
end

function RTs_SE(a,b,distri,cohmeasure)
    # a quantity to mean
    # b quantity to test
    # Error after Sucess
    t = 0
    c = 0
    output=Float64[]
    for i=1:length(a)-1
        x = b[i]
        y=b[i+1]
        x > 0 && y <0 && distri[i+1]*distri[i] > 0 && distri[i+1] == cohmeasure || continue
        push!(output,a[i+1])

    end
    return mean(output),std(output),length(output)
end

function RTs_ES(a,b,distri,cohmeasure)
    # a quantity to mean
    # b quantity to test
    ## Sucess after Error
    output=Float64[]
    for i=1:length(a)-1
        x = b[i]
        y=b[i+1]
        x < 0 && y >0 && distri[i+1] == cohmeasure|| continue

        push!(output,a[i+1])

    end
    return mean(output),std(output),length(output)
end

function RTs_SSS(a,b)
    # Sucess after Sucess after Sucess
    # same sign coh but not same one
    output=Float64[]
    for i=1:length(a)-2
        x = b[i]
        y=b[i+1]
        z=b[i+2]
        x > 0 && y >0 && z>0 && distri[i]*distri[i+1]>0 && distri[i+1]*distri[i+2] && distri[i+2] == cohmeasure|| continue
        push!(output,a[i+2])
    end
    return mean(output),std(output),length(output)
end

function RTs_SSE(a,b)
    # Error after Sucess after Sucess
    t = 0
    c = 0
    output=Float64[]

    for i=1:length(a)-2
        x = b[i]
        y=b[i+1]
        z=b[i+2]
        x > 0 && y >0 && z<0&& distri[i]*distri[i+1]>0 && distri[i+1]*distri[i+2] && distri[i+2] == cohmeasure|| continue
        push!(output,a[i+2])
    end
    return mean(output),std(output),length(output)
end

function RTs_SS_chgt(RT,Res,distri,coherence)
    output = Float64[]
    for i =1:length(RT)-1
        if Res[i]==1 && Res[i+1] == 1 && distri[i+1] == coherence && distri[i]*distri[i+1] <0  || continue
            push!(output,RT[i+1])
        end
    end
    return mean(output)
end

function analysis_RTs_E(result,Rts_coh)
        listc=Float64[]
        listv=Float64[]

        for c in [0 1 2 3 4 5 6 7 8]
                tempv = Float64[]
                temp=Float64[]
                for i in 1:length(result["distri"])-1
                        if sign(result["distri"][i])==sign(result["distri"][i+1]) &&
                                result["Confidence"][i+1]>(c-1) &&result["NotError"][i]==-1 #&&result["NotError"][i+1]==1
                                x=result["distri"][i+1]
                                push!(tempv,result["RTs"][i+1])
                                if abs(x) <10
                                push!(temp,(Rts_coh["$x"]-result["RTs"][i+1])/Rts_coh["$x"])
end
                        end
                end
                push!(listc,mean(temp))
                push!(listv,var(temp))
        end
        return listc,listv
end

function analysis_RTs_SE(result,Rts_coh)
    listc=Float64[]
    listv=Float64[]

    for c in [0 1 2 3 4 5 6 7 8]
            tempv = Float64[]
            temp=Float64[]
            for i in 1:length(result["distri"])-1
                    if sign(result["distri"][i])==sign(result["distri"][i+1]) &&
                            result["Confidence"][i+1]>(c-1) &&result["NotError"][i]==-1 &&result["NotError"][i+1]==1
                            x=result["distri"][i+1]
                            push!(tempv,result["RTs"][i+1])
                            push!(temp,(Rts_coh["$x"]-result["RTs"][i+1])/Rts_coh["$x"])
                    end
            end
            push!(listc,mean(temp))
            push!(listv,var(temp))
    end
    return listc,listv
end

function analysis_RTs_SS(result,Rts_coh)
    listc=Float64[]
    listv=Float64[]

    for c in [0 1 2 3 4 5 6 7 8]
            tempv = Float64[]
            temp=Float64[]
            for i in 1:length(result["distri"])-1
                    if sign(result["distri"][i])==sign(result["distri"][i+1]) &&
                            result["Confidence"][i+1]>(c-1) &&result["NotError"][i]==1 &&result["NotError"][i+1]==1
                            x=result["distri"][i+1]
                            push!(tempv,result["RTs"][i+1])
                            push!(temp,(Rts_coh["$x"]-result["RTs"][i+1])/Rts_coh["$x"])
                    end
            end
            push!(listc,mean(temp))
            push!(listv,var(temp))
    end
    return listc,listv
end

function analysis_RTs_E_cohspec(result,Rts_coh)
        listc=Float64[]
        for c in [0 1 2 3 4 5 6 7 8]
                temp=Float64[]
                for i in 1:length(result["distri"])-1
                        if sign(result["distri"][i])==sign(result["distri"][i+1]) &&
                                result["Confidence"][i+1]==c &&result["NotError"][i]==-1 #&&result["NotError"][i+1]==1
                                x=result["distri"][i+1]
                                push!(temp,(Rts_coh["$x"]-result["RTs"][i+1])/Rts_coh["$x"])
                        end
                end
                push!(listc,mean(temp))
        end
        return listc
end

function analysis_RTs_SE_cohspec(result,Rts_coh)
        listc=Float64[]
        for c in [0 1 2 3 4 5 6 7 8]
                temp=Float64[]
                for i in 1:length(result["distri"])-1
                        if sign(result["distri"][i])==sign(result["distri"][i+1]) &&
                                result["Confidence"][i+1]==c &&result["NotError"][i]==-1 &&result["NotError"][i+1]==1
                                x=result["distri"][i+1]
                                push!(temp,(Rts_coh["$x"]-result["RTs"][i+1])/Rts_coh["$x"])
                        end
                end
                push!(listc,mean(temp))
        end
        return listc
end

function analysis_RTs_SS_cohspec(result,Rts_coh)
        listc=Float64[]
        for c in [0 1 2 3 4 5 6 7 8]
                temp=Float64[]
                for i in 1:length(result["distri"])-1
                        if sign(result["distri"][i])==sign(result["distri"][i+1]) &&
                                result["Confidence"][i+1]==(c-1) &&result["NotError"][i]==1 &&result["NotError"][i+1]==1
                                x=result["distri"][i+1]
                                push!(temp,(Rts_coh["$x"]-result["RTs"][i+1])/Rts_coh["$x"])
                        end
                end
                push!(listc,mean(temp))
        end
        return listc
end

function load_Jerome_E()

       data = CSV.read("/home/kevin/Documents/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)
        #data = CSV.read("/media/kevin/Boulot/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)

        rt_confidence = [] # for confidence task
        results_confidence = []
        cursor_confidence = []
        coherence_confidence = []
        reponse_confidence = []
        name_confidence=[]
        number=Float64[]

        for i=1:length(data[1])
                if  data[3][i] == "confidence"
                        push!(reponse_confidence,data[6][i])
                        push!(cursor_confidence,data[11][i])
                        push!(coherence_confidence,data[5][i])
                        push!(rt_confidence,data[7][i])
                        push!(results_confidence,data[20][i])
                        push!(name_confidence,data[2][i])
                        push!(number,parse(Float64,data[1][i]))
                    end
        end
        #list_coherence = [-2.4 -1.6 -0.8 0 0.8 1.6 2.4]
        #list_coherence_string=["-2.4" "-1.6" "-0.8" "0" "0.8" "1.6" "2.4"]
        list_coherence = [-2.4 -1.6 -0.8 0 0.8 1.6 2.4]
        list_coherence_string=["-2.4" "-1.6" "-0.8" "0" "0.8" "1.6" "2.4" ]

        dic_rt_confidence = Dict()
        for coh_string in list_coherence_string
                rt_new = rt_confidence[coherence_confidence.==parse(Float64,coh_string)]
                dic_rt_confidence[coh_string]=mean(rt_new[rt_new.<1500])
        end
        dic_rt_cursor_error=Dict()
        for cursor in [0 1 2 3 4 5 6 7 8]
             rt_list_difference=Float64[]

             for i=1:(length(rt_confidence)-1)
                if  sign(coherence_confidence[i])==sign(coherence_confidence[i+1])  &&
                    name_confidence[i] == name_confidence[i+1] &&
                    results_confidence[i] == 0 && number[i] +1 == number[i+1]#&&results_confidence[i+1] == 1

                   if rt_confidence[i]<1500 && rt_confidence[i+1] <1500 &&
                      cursor_confidence[i]>cursor #&& cohVEA[i]!=0

                        x=coherence_confidence[i+1]
                        if x==0 # because 0 is a float in our data
                            push!(rt_list_difference,((dic_rt_confidence["0"]-rt_confidence[i+1])/dic_rt_confidence["0"]))
                        else
                            push!(rt_list_difference,((dic_rt_confidence["$x"]-rt_confidence[i+1])/dic_rt_confidence["$x"]))

                    end
                    end
                end
             end
    println(length(rt_list_difference))
    dic_rt_cursor_error["$cursor"]=mean(rt_list_difference)*100
    dic_rt_cursor_error["var $cursor"]=var(rt_list_difference*100)/700
    dic_rt_cursor_error["label $cursor"]=length(rt_list_difference)
    end

    exp_list2=Float64[]
    var_list2=Float64[]
    for cursor in [0 1 2 3 4 5 6 7 8]
                push!(exp_list2,dic_rt_cursor_error["$cursor"])
                push!(var_list2,dic_rt_cursor_error["var $cursor"])
    end
    return exp_list2,var_list2
end

function load_Jerome_SS()

            data = CSV.read("/home/kevin/Documents/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)
        #    data = CSV.read("/media/kevin/Boulot/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)
            rt_confidence = [] # for confidence task
        results_confidence = []
        cursor_confidence = []
        coherence_confidence = []
        reponse_confidence = []
        name_confidence=[]
        number=[]
        for i=1:length(data[1])
                if  data[3][i] == "confidence"
                        push!(reponse_confidence,data[6][i])
                        push!(cursor_confidence,data[11][i])
                        push!(coherence_confidence,data[5][i])
                        push!(rt_confidence,data[7][i])
                        push!(results_confidence,data[20][i])
                        push!(name_confidence,data[2][i])
                        push!(number,parse(Float64,data[1][i]))
                    end
        end
        list_coherence = [-2.4 -1.6 -0.8 0 0.8 1.6 2.4]
        list_coherence_string=["-2.4" "-1.6" "-0.8" "0" "0.8" "1.6" "2.4"]
        dic_rt_confidence = Dict()
        for coh_string in list_coherence_string
                rt_new = rt_confidence[coherence_confidence.==parse(Float64,coh_string)]
                dic_rt_confidence[coh_string]=mean(rt_new[rt_new.<1500])
                println(mean(rt_new[rt_new.<1500]))
        end
        dic_rt_cursor_error=Dict()
        for cursor in [0 1 2 3 4 5 6 7 8]
             rt_list_difference=Float64[]

             for i=1:(length(rt_confidence)-1)
                if  coherence_confidence[i]==coherence_confidence[i+1]  &&
                    name_confidence[i] == name_confidence[i+1] &&
                    number[i]+1==number[i+1] &&
                    results_confidence[i] == 1 &&results_confidence[i+1] == 1

                   if rt_confidence[i]<1500 && rt_confidence[i+1] <1500 &&
                      cursor_confidence[i]>cursor #&& cohVEA[i]!=0

                        x=coherence_confidence[i+1]
                        if x==0 # because 0 is a float in our data
                            push!(rt_list_difference,((rt_confidence[i]-rt_confidence[i+1])/dic_rt_confidence["0"]))
                        else
                            push!(rt_list_difference,((rt_confidence[i]-rt_confidence[i+1])/dic_rt_confidence["$x"]))
                        end

                        #if x==0 # because 0 is a float in our data
                        #    push!(rt_list_difference,((dic_rt_confidence["0"]-rt_confidence[i+1])/dic_rt_confidence["0"]))
                        #else
                        #    push!(rt_list_difference,((dic_rt_confidence["$x"]-rt_confidence[i+1])/dic_rt_confidence["$x"]))

                    #end
                    end
                end
             end
    println(length(rt_list_difference))
    dic_rt_cursor_error["$cursor"]=mean(rt_list_difference)*100
    dic_rt_cursor_error["var $cursor"]=var(rt_list_difference)
    dic_rt_cursor_error["label $cursor"]=length(rt_list_difference)
    end

    exp_list2=Float64[]
    var_list2=Float64[]
    for cursor in [0 1 2 3 4 5 6 7 8]
                push!(exp_list2,dic_rt_cursor_error["$cursor"])
                push!(var_list2,dic_rt_cursor_error["var $cursor"])
    end
    return exp_list2,var_list2
end

function load_Jerome_SE_cohspec()

            data = CSV.read("/home/kevin/Documents/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)
        #    data = CSV.read("/media/kevin/Boulot/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)
    rt_confidence = [] # for confidence task
        results_confidence = []
        cursor_confidence = []
        coherence_confidence = []
        reponse_confidence = []
        name_confidence=[]

        for i=1:length(data[1])
                if  data[3][i] == "confidence"
                        push!(reponse_confidence,data[6][i])
                        push!(cursor_confidence,data[11][i])
                        push!(coherence_confidence,data[5][i])
                        push!(rt_confidence,data[7][i])
                        push!(results_confidence,data[20][i])
                        push!(name_confidence,data[2][i])
                    end
        end
        list_coherence = [-2.4 -1.6 -0.8 0 0.8 1.6 2.4]
        list_coherence_string=["-2.4" "-1.6" "-0.8" "0" "0.8" "1.6" "2.4"]
        dic_rt_confidence = Dict()
        for coh_string in list_coherence_string
                rt_new = rt_confidence[coherence_confidence.==parse(Float64,coh_string)]
                dic_rt_confidence[coh_string]=mean(rt_new[rt_new.<1500])
        end
        dic_rt_cursor_error=Dict()
        for cursor in [0 1 2 3 4 5 6 7 8]
             rt_list_difference=Float64[]

             for i=1:(length(rt_confidence)-1)
                if  sign(coherence_confidence[i])==sign(coherence_confidence[i+1])  &&
                    name_confidence[i] == name_confidence[i+1] &&
                    results_confidence[i] == 0 &&results_confidence[i+1] == 1

                   if rt_confidence[i]<1500 && rt_confidence[i+1] <1500 &&
                      cursor_confidence[i]==cursor+1 #&& cohVEA[i]!=0

                        x=coherence_confidence[i+1]
                        if x==0 # because 0 is a float in our data
                            push!(rt_list_difference,((dic_rt_confidence["0"]-rt_confidence[i+1])/dic_rt_confidence["0"]))
                        else
                            push!(rt_list_difference,((dic_rt_confidence["$x"]-rt_confidence[i+1])/dic_rt_confidence["$x"]))
                        end
                    end
                end
             end
    println(length(rt_list_difference))
    dic_rt_cursor_error["$cursor"]=mean(rt_list_difference)*100
    dic_rt_cursor_error["var $cursor"]=var(rt_list_difference)
    dic_rt_cursor_error["label $cursor"]=length(rt_list_difference)
    end

    exp_list2=Float64[]
    var_list2=Float64[]
    for cursor in [0 1 2 3 4 5 6 7 8]
                push!(exp_list2,dic_rt_cursor_error["$cursor"])
                push!(var_list2,dic_rt_cursor_error["var $cursor"])
    end
    return exp_list2,var_list2
end

function load_Jerome_SE()

            data = CSV.read("/home/kevin/Documents/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)
        #    data = CSV.read("/media/kevin/Boulot/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)
rt_confidence = [] # for confidence task
        results_confidence = []
        cursor_confidence = []
        coherence_confidence = []
        reponse_confidence = []
        name_confidence=[]
        number=Float64[]
        for i=1:length(data[1])
                if  data[3][i] == "confidence"
                        push!(reponse_confidence,data[6][i])
                        push!(cursor_confidence,data[11][i])
                        push!(coherence_confidence,data[5][i])
                        push!(rt_confidence,data[7][i])
                        push!(results_confidence,data[20][i])
                        push!(name_confidence,data[2][i])
                        push!(number,parse(Float64,data[1][i]))
                    end
        end
        list_coherence = [-2.4 -1.6 -0.8 0 0.8 1.6 2.4]
        list_coherence_string=["-2.4" "-1.6" "-0.8" "0" "0.8" "1.6" "2.4"]
        dic_rt_confidence = Dict()
        for coh_string in list_coherence_string
                rt_new = rt_confidence[coherence_confidence.==parse(Float64,coh_string)]
                dic_rt_confidence[coh_string]=mean(rt_new[rt_new.<1500])
        end
        dic_rt_cursor_error=Dict()
        for cursor in [0 1 2 3 4 5 6 7 8]
             rt_list_difference=Float64[]

             for i=1:(length(rt_confidence)-1)
                if  sign(coherence_confidence[i])==sign(coherence_confidence[i+1])  &&
                    name_confidence[i] == name_confidence[i+1] &&
                    number[i] + 1 == number[i+1] &&
                    results_confidence[i] == 0 &&results_confidence[i+1] == 1

                   if rt_confidence[i]<1500 && rt_confidence[i+1] <1500 &&
                      cursor_confidence[i]>=cursor+1 #&& cohVEA[i]!=0

                        x=coherence_confidence[i+1]
                        if x==0 # because 0 is a float in our data
                            push!(rt_list_difference,((dic_rt_confidence["0"]-rt_confidence[i+1])/dic_rt_confidence["0"]))
                        else
                            push!(rt_list_difference,((dic_rt_confidence["$x"]-rt_confidence[i+1])/dic_rt_confidence["$x"]))
                        end
                    end
                end
             end
    println(length(rt_list_difference))
    dic_rt_cursor_error["$cursor"]=mean(rt_list_difference)*100
    dic_rt_cursor_error["var $cursor"]=var(rt_list_difference)
    dic_rt_cursor_error["label $cursor"]=length(rt_list_difference)
    end

    exp_list2=Float64[]
    var_list2=Float64[]
    for cursor in [0 1 2 3 4 5 6 7 8]
                push!(exp_list2,dic_rt_cursor_error["$cursor"])
                push!(var_list2,dic_rt_cursor_error["var $cursor"])
    end
    return exp_list2,var_list2
end

function load_Jerome_nextRTE()

            data = CSV.read("/home/kevin/Documents/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)
        #    data = CSV.read("/media/kevin/Boulot/PhD/PhD_Data/Jerome_17/secondFirst.csv",nullable=false)
rt_confidence = [] # for confidence task
        results_confidence = []
        cursor_confidence = []
        coherence_confidence = []
        reponse_confidence = []
        name_confidence=[]
        number=Float64[]
        for i=1:length(data[1])
                if  data[3][i] == "confidence"
                        push!(reponse_confidence,data[6][i])
                        push!(cursor_confidence,data[11][i])
                        push!(coherence_confidence,data[5][i])
                        push!(rt_confidence,data[7][i])
                        push!(results_confidence,data[20][i])
                        push!(name_confidence,data[2][i])
                        push!(number,parse(Float64,data[1][i]))
                    end
        end
        list_coherence = [-2.4 -1.6 -0.8 0 0.8 1.6 2.4]
        list_coherence_string=["-2.4" "-1.6" "-0.8" "0" "0.8" "1.6" "2.4"]
        dic_rt_confidence = Dict()
        for coh_string in list_coherence_string
                rt_new = rt_confidence[coherence_confidence.==parse(Float64,coh_string)]
                dic_rt_confidence[coh_string]=mean(rt_new[rt_new.<1500])
        end
        dic_rt_cursor_error=Dict()
        for cursor in [0 1 2 3 4 5 6 7 8]
             rt_list_difference=Float64[]

             for i=1:(length(rt_confidence)-1)
                if  name_confidence[i] == name_confidence[i+1] &&
                    number[i] + 1 == number[i+1] &&
                    results_confidence[i] == 0 # &&results_confidence[i+1] == 1

                   if rt_confidence[i]<1500 && rt_confidence[i+1] <1500 &&
                      cursor_confidence[i]==cursor+1 #&& cohVEA[i]!=0

                        x=coherence_confidence[i+1]
                        if x==0 # because 0 is a float in our data
                            push!(rt_list_difference,((dic_rt_confidence["0"]-rt_confidence[i+1])))
                        else
                            push!(rt_list_difference,((dic_rt_confidence["$x"]-rt_confidence[i+1])))
                        end
                    end
                end
             end
    println(length(rt_list_difference))
    dic_rt_cursor_error["$cursor"]=mean(rt_list_difference)
    dic_rt_cursor_error["var $cursor"]=var(rt_list_difference)
    dic_rt_cursor_error["label $cursor"]=length(rt_list_difference)
    end

    exp_list2=Float64[]
    var_list2=Float64[]
    for cursor in [0 1 2 3 4 5 6 7 8]
                push!(exp_list2,dic_rt_cursor_error["$cursor"])
                push!(var_list2,dic_rt_cursor_error["var $cursor"])
    end
    return exp_list2,var_list2
end


text:
module BartIO

using BufferedStreams
using ConfParser
using PyCall

# Exported functions
export readcfl
export writecfl
export initBart

"""
    bart = initBart(pathtobart::String,pathtobartpy::String)

Initialize the installation of bart and to bartpy in order to make it available
from the bartpy package through PyCall and store the path in a config file

## Input Parameters
-

## output

# Example
"""

function initBart(path2bart::String="",path2bartpy::String="")

    conf = ConfParser.ConfParse("confs/config.ini")
    parse_conf!(conf)

    pathtobart=CheckAndSetPath!(conf,"BART","pathtobart",path2bart)
    pathtobartpy=CheckAndSetPath!(conf,"BART","pathtobartpy",path2bartpy)

    # Build PyBart

    path2BartPython = pathtobart*"/python"
    py"""
    import sys
    import os
    sys.path.insert(0, $path2BartPython)
    os.environ['TOOLBOX_PATH'] = $pathtobart
    """
    python_pycall = PyCall.python

    cmd = `cd $pathtobartpy \; $python_pycall setup.py install`
    run(cmd)

    #@PyCall.pyimport bartpy.tools as bartpy #Equivalent to -> bartpy = pyimport("bartpy.tools") but does not work in module...
    bartpy = pyimport("bartpy.tools")
    bartpy.version()

    return bartpy
end

## Utility functions
function CheckAndSetPath!(conf::ConfParse,blockname::String,pathname::String,path::String)

    if isempty(path);
        path = retrieve(conf,blockname,pathname)
        if isempty(path); error("$pathname is not defined ! Set it with the function : \n initBart(path2bart::String,path2bartpy::String)"); end
    else
        commit!(conf, blockname, pathname, path);
        save!(conf)
    end

    # check if path exists
    if !isdir(path)
        error(path*" is not a valid directory ! redefined it")
    end
    return path
end

"""
    readcfl(filename::String)

- readcfl(filename(no extension)) -> Array{ComplexF32,N} where N is defined the filename.hdr file
- readcfl(filename.cfl) -> Array{ComplexF32,N} where N is defined the filename.hdr file
- readcfl(filename.hdr) -> Array{ComplexF32,N} where N is defined the filename.hdr file

Reads complex data from files created by the Berkeley Advanced Reconstruction Toolbox (BART).
The output is an Array of ComplexF32 with the dimensions stored in a .hdr file.

## Parameters:
- filename:   path and filename of the cfl and hdr files, which can either be without extension, end on .cfl, or end on .hdr
"""
function readcfl(filename::String)

    if filename[end-3:end] == ".cfl"
        filenameBase = filename[1:end-4]
    elseif filename[end-3:end] == ".hdr"
        filenameBase = filename[1:end-4]
        filename = string(filenameBase, ".cfl");
    else
        filenameBase = filename
        filename = string(filenameBase, ".cfl");
    end

    dims = readreconheader(filenameBase);
    data = Array{ComplexF32}(undef, Tuple(dims))

    fid = BufferedInputStream(open(filename));

    for i in eachindex(data)
        data[i] = read(fid, Float32) + 1im * read(fid, Float32)
    end

    close(fid);
    return data
end

function readreconheader(filenameBase::String)
    filename = string(filenameBase, ".hdr");
    fid = open(filename);

    line = ["#"]
    while line[1] == "#"
        line = split(readline(fid))
    end

    dims = parse.(Int, line)
    close(fid);
    return dims
end

"""
    writecfl(filename::String,dataCfl::Array{ComplexF32})

- writecfl(filename(no extension),Array{ComplexF32})
- writecfl(filename.cfl, Array{ComplexF32})
- writecfl(filename.hdr,Array{ComplexF32})

Write complex data to files following the convention of the Berkeley Advanced Reconstruction Toolbox (BART).
The input is an Array of ComplexF32 with the dimensions stored in a .hdr file.

## Parameters:
- filename:   path and filename of the cfl and hdr files, which can either be without extension, end on .cfl, or end on .hdr
- Array{ComplexF32,N}:   Array of ComplexF32 corresponding to image/k-space

"""
function writecfl(filename::String,dataCfl::Array{ComplexF32})

    if filename[end-3:end] == ".cfl"
        filenameBase = filename[1:end-4]
    elseif filename[end-3:end] == ".hdr"
        filenameBase = filename[1:end-4]
        filename = string(filenameBase, ".cfl");
    else
        filenameBase = filename
        filename = string(filenameBase, ".cfl");
    end

    dimTuple = size(dataCfl)
    dims = ones(Int,16,1)

    for i in 1:length(dimTuple)
        dims[i]=dimTuple[i];
    end

    writereconheader(filenameBase,dims);

    fid = BufferedOutputStream(open(filename,"w"));
    write(fid,dataCfl)
    close(fid);
end

function writereconheader(filenameBase::String,dims::Array{Int})
    filename = string(filenameBase, ".hdr");

    fid = open(filename,"w");
    write(fid,"# Dimension\n")
    a = length(dims)
    for i in 1:length(dims)
        write(fid,string(dims[i])*" ")
    end
    close(fid)
end

end # module

ext: jl


text:
module BenchTuple
using BenchmarkTools
using Statistics
using Gadfly
using MLStyle
using DataFrames
import Match
import Rematch
using ..ArbitrarySampler
using ..Utils

mod′(n) = x -> mod(x, n)
len_in(rng) = x -> length(x) in rng
is_in(rng) = x -> x in rng

integer20 = @spec ::Integer isa mod′(20)

struct UserTy{T}
    a::T
    b::Symbol
end
@as_record UserTy

specs = [
    :spec1 => (@spec (1, (_, "2", _), ("3", 4, 5))),
    :spec2 => (@spec (_, "1", 2, _, (3, "4", _), _)),
    :spec3 => (@spec (_, 1, _, 2, _, 3, _, 4, _, 5)),
    :spec4 => (@spec ((1, 2, 3, _), (4, 5, 6, _, (7, 8, 9, _, (11, 12, 13))))),
    :spec5 => (@spec (::String, ::Symbol, ::Real, ([1, _], UserTy(:a, :b)))),
    :_ => @spec(_),
]

implementations = [
    :MLStyle => (@λ begin
        (1, (_, "2", _), ("3", 4, 5)) -> 1
        (_, "1", 2, _, (3, "4", _), _) -> 2
        (_, 1, _, 2, _, 3, _, 4, _, 5) -> 3
        ((1, 2, 3, _), (4, 5, 6, _, (7, 8, 9, _, (11, 12, 13)))) -> 4
        _ -> 5
    end),
    :Rematch => function (x)
        Rematch.@match x begin
            (1, (_, "2", _), ("3", 4, 5)) => 1
            (_, "1", 2, _, (3, "4", _), _) => 2
            (_, 1, _, 2, _, 3, _, 4, _, 5) => 3
            ((1, 2, 3, _), (4, 5, 6, _, (7, 8, 9, _, (11, 12, 13)))) => 4
            _ => 5
        end
    end,
    Symbol("Match.jl") => function (x)
        Match.@match x begin
            (1, (_, "2", _), ("3", 4, 5)) => 1
            (_, "1", 2, _, (3, "4", _), _) => 2
            (_, 1, _, 2, _, 3, _, 4, _, 5) => 3
            ((1, 2, 3, _), (4, 5, 6, _, (7, 8, 9, _, (11, 12, 13)))) => 4
            _ => 5
        end
    end,
    :HandWritten => function (tp)
        !(tp isa Tuple) ? 5 :
                let n = length(tp)
            if n === 3 && tp[1] === 1
                tp_ = tp[2]
                if !(tp isa Tuple) || tp_[2] != "2"
                    return 5
                end
                tp[3] == ("3", 4, 5) ? 1 : 5
            elseif n === 6 && tp[2] == "1" && tp[3] === 2
                tp = tp[5]
                !(tp isa Tuple) ? 5 : tp[1] === 3 && tp[2] == "4" ? 2 : 5
            elseif n === 10 &&
                   tp[2] === 1 &&
                   tp[4] === 2 &&
                   tp[6] === 3 &&
                   tp[8] === 4 &&
                   tp[10] === 5
                3
            elseif n === 2 && tp[1] isa Tuple && tp[2] isa Tuple
                @inline eqtp(a, slice, v) =
                    all(slice) do i
                        a[i] === v[i]
                    end
                (a, b) = tp
                if eqtp(a, 1:3, (1, 2, 3)) &&
                   eqtp(b, 1:3, (4, 5, 6)) &&
                   let a = b[5]
                       eqtp(a, 1:3, (7, 8, 9)) && a[5] === (11, 12, 13)
                   end
                    4
                end
            else
                5
            end
        end
    end,
]

records = NamedTuple{(:time_mean, :implementation, :case)}[]
for (spec_id, spec) in specs
    # group_key = string(spec_id)
    # suite[group_key] = BenchmarkGroup()
    for (impl_id, impl_fn) in implementations
        bench′ =
            @benchmark $impl_fn(sample) setup = (sample = $generate($spec)) samples = 2000
        time′ = mean(bench′.times)
        @info :bench (spec_id, impl_id, time′)
        push!(records, (time_mean = time′, implementation = impl_id, case = spec_id))
    end
end

df = DataFrame(records)

@info df

theme = Theme(
    guide_title_position = :left,
    colorkey_swatch_shape = :circle,
    minor_label_font = "Consolas",
    major_label_font = "Consolas",
    point_size = 5px,
)

report_meantime, df_time = report(df, Guide.title("Tuples"); benchfield = :time_mean)

open("stats/bench-tuple.txt", "w") do f
    write(f, string(df))
end

draw(SVG("stats/bench-tuple.svg", 14inch, 6inch), report_meantime)

end


text:
module TestClassifier
	using Test
	using kNN
	using DataFrames
    using RDatasets
    using Distances
    using StatsBase
	using Statistics

    iris = dataset("datasets", "iris")
    X = collect(convert(Array, iris[1:4])')
    y = convert(Array, iris[5])
    model = knn(X, y, metric = Euclidean())

    predict_k1 = predict(model, X, 1)
    predict_k2 = predict(model, X, 2)
    predict_k3 = predict(model, X, 3)
    predict_k4 = predict(model, X, 4)
    predict_k5 = predict(model, X, 5)

    mean(predict_k1 .== y)
    mean(predict_k2 .== y)
    mean(predict_k3 .== y)
    mean(predict_k4 .== y)
    mean(predict_k5 .== y)
end


text:
module Semagrams

include("SVG.jl")
include("Boxes.jl")
include("Logo.jl")
include("Muesli.jl")
include("Schema.jl")
include("Data.jl")
include("UI.jl")
include("Examples.jl")

using Reexport

@reexport using .SVG
@reexport using .Boxes
@reexport using .Logo
@reexport using .Muesli
@reexport using .Schema
@reexport using .Data
@reexport using .UI
# Examples are *not* reexported by default,
# because they could clash with similar definitions elsewhere

end


text:
module Pkg3Tests

include("pkg.jl")
include("resolve.jl")

end # module

max_stars_repo_name: djsegal/julia-fork

max_issues_repo_name: djsegal/julia-fork

max_forks_repo_name: djsegal/julia-fork


text:
struct ICKey{IC <: InitialConditionType, D <: PSY.Component}
    ic_type::Type{IC}
    device_type::Type{D}
end

struct InitialConditions
    use_parameters::Bool
    data::Dict{ICKey, Vector{InitialCondition}}
end

function InitialConditions(;
    use_parameters = false,
    data = Dict{ICKey, Array{InitialCondition}}(),
)
    return InitialConditions(use_parameters, data)
end

get_use_parameters(container::InitialConditions) = container.use_parameters

function has_initial_conditions(container::InitialConditions, key::ICKey)
    return key in keys(container.data)
end

function get_initial_conditions(container::InitialConditions, key::ICKey)
    initial_conditions = get(container.data, key, nothing)
    if initial_conditions === nothing
        throw(IS.InvalidValue("initial conditions are not stored for $(key)"))
    end

    return initial_conditions
end

function set_initial_conditions!(container::InitialConditions, key::ICKey, value)
    @debug "set_initial_condition_container" key
    container.data[key] = value
end

"""
Iterate over the keys and vectors of initial conditions.
"""
iterate_initial_conditions(container::InitialConditions) = pairs(container.data)


text:
function Epoch{S}(date::Date, time::Time, args...) where S
    daysec = round(Int64, (j2000(date) - 0.5) * SECONDS_PER_DAY)
    horsec = round(Int64, hour(time) * SECONDS_PER_HOUR)
    minsec = round(Int64, minute(time) * SECONDS_PER_MINUTE)
    sec = daysec + horsec + minsec + second(time)

    return Epoch{S}(sec, time.fraction)
end

Dates.default_format(::Type{<:Epoch}) = EPOCH_ISO_FORMAT[]

"""
    Epoch(str[, format])

Construct an `Epoch` from a string `str`. Optionally a `format` definition can
be passed as a [`DateFormat`](https://docs.julialang.org/en/stable/stdlib/Dates/#Dates.DateFormat)
object or as a string. In addition to the character codes supported by `DateFormat` the character
code `D` is supported which is parsed as "day of year" (see the example below) and the character
code `t` which is parsed as the time scale.  The default format is `yyyy-mm-ddTHH:MM:SS.sss ttt`.

**Note:** Please be aware that this constructor requires that the time scale is part of `str`, e.g.
`2018-02-06T00:00 TAI`. Otherwise use an explicit constructor, e.g. `Epoch{TAI}`.

### Example ###

```jldoctest; setup = :(using AstroTime)
julia> Epoch("2018-02-06T20:45:00.0 TAI")
2018-02-06T20:45:00.000 TAI

julia> Epoch("2018-037T00:00 TAI", "yyyy-DDDTHH:MM ttt")
2018-02-06T00:00:00.000 TAI
```
"""
Epoch(str::AbstractString, format::Dates.DateFormat=Dates.default_format(Epoch)) = parse(Epoch, str, format)

Epoch(str::AbstractString, format::AbstractString) = Epoch(str, Dates.DateFormat(format))

"""
    Epoch{S}(str[, format]) where S

Construct an `Epoch` with time scale `S` from a string `str`. Optionally a `format` definition can
be passed as a [`DateFormat`](https://docs.julialang.org/en/stable/stdlib/Dates/#Dates.DateFormat)
object or as a string. In addition to the character codes supported by `DateFormat` the code `D` can
be used which is parsed as "day of year" (see the example below).  The default format is
`yyyy-mm-ddTHH:MM:SS.sss`.

### Example ###

```jldoctest; setup = :(using AstroTime)
julia> Epoch{InternationalAtomicTime}("2018-02-06T20:45:00.0")
2018-02-06T20:45:00.000 TAI

julia> Epoch{InternationalAtomicTime}("February 6, 2018", "U d, y")
2018-02-06T00:00:00.000 TAI

julia> Epoch{InternationalAtomicTime}("2018-037T00:00", "yyyy-DDDTHH:MM")
2018-02-06T00:00:00.000 TAI
```
"""
function Epoch{S}(str::AbstractString, format::Dates.DateFormat=Dates.default_format(Epoch{S})) where S
    return Dates.parse(Epoch{S}, str, format)
end

Epoch{S}(str::AbstractString, format::AbstractString) where {S} = Epoch{S}(str, Dates.DateFormat(format))

Epoch{S}(d::Date) where {S} = Epoch{S}(d, AstroDates.H00)

Epoch{S}(dt::DateTime) where {S} = Epoch{S}(Date(dt), Time(dt))
Epoch{S}(dt::Dates.DateTime) where {S} = Epoch{S}(DateTime(dt))

"""
    Epoch{S}(year, month, day, hour=0, minute=0, second=0.0) where S

Construct an `Epoch` with time scale `S` from date and time components.

### Example ###

```jldoctest; setup = :(using AstroTime)
julia> Epoch{InternationalAtomicTime}(2018, 2, 6, 20, 45, 0.0)
2018-02-06T20:45:00.000 TAI

julia> Epoch{InternationalAtomicTime}(2018, 2, 6)
2018-02-06T00:00:00.000 TAI
```
"""
function Epoch{S}(year::Int, month::Int, day::Int, hour::Int=0,
                  minute::Int=0, second::Float64=0.0, args...) where S
    Epoch{S}(Date(year, month, day), Time(hour, minute, second), args...)
end

function Epoch{S}(year::Int64, month::Int64, day::Int64, dayofyear::Int64,
                  hour::Int64, minute::Int64, second::Int64, milliseconds::Int64,
                  fractionofsecond::Float64) where S
    dt = AstroDates.DateTime(year, month, day, dayofyear, hour, minute, second,
                             milliseconds, fractionofsecond)
    return Epoch{S}(dt)
end

function Epoch(year::Int64, month::Int64, day::Int64, dayofyear::Int64,
               hour::Int64, minute::Int64, second::Int64, milliseconds::Int64,
               fractionofsecond::Float64, scale::S) where S<:TimeScale
    if scale === TimeScales.NotATimeScale()
        throw(ArgumentError("Could not parse the provided string as an `Epoch`." *
                            " No time scale was provided."))
    end

    return Epoch{S}(year, month, day, dayofyear,
                    hour, minute, second, milliseconds,
                    fractionofsecond)
end

max_stars_repo_name: JuliaAstrodynamics/AstronomicalTime

max_issues_repo_name: JuliaAstrodynamics/AstronomicalTime

max_forks_repo_name: JuliaAstrodynamics/AstronomicalTime


text:
function inv(A)
    hasInv = issingular_convert!(copy(A))

    if hasInv == false
        return nothing
    end

    rows, columns = size(A)
    E = A^0
    D = Matrix{eltype(A)}(undef, rows, 2columns)

    for i in 1:rows
        D[i, :] = append!(copy(A[i, :]), E[i, :])
    end

    convert!(D)

    for i in 1:rows
        E[i, :] = D[i, rows+1:2rows]
    end

    return E
end


text:
# This file is a part of ValueShapes.jl, licensed under the MIT License (MIT).

using ValueShapes
using Test

using Random
using ElasticArrays
using ArraysOfArrays
using FillArrays
using InverseFunctions, ChangesOfVariables
using ChainRulesCore: rrule, NoTangent
import TypedTables

@testset "abstract_value_shape" begin
    @testset "default_datatype" begin
        @test @inferred(ValueShapes.default_datatype(Integer)) == Int
        @test @inferred(ValueShapes.default_datatype(AbstractFloat)) == Float64
        @test @inferred(ValueShapes.default_datatype(Real)) == Float64
        @test @inferred(ValueShapes.default_datatype(Complex)) == Complex{Float64}

        @test @inferred(elshape(Complex(1.0, 2.0))) == ScalarShape{Complex{Float64}}()
        @test @inferred(elshape([[3, 5], [3, 2]])) == ArrayShape{Int,1}((2,))

        @test ValueShapes.stripscalar(Ref(ScalarShape{Real})) == ScalarShape{Real}

        # @test Vector{Real}(undef, ArrayShape{Real}((2,1))) ==  # weird typing going on with default_unshaped_eltype

        arrshape = ArrayShape{Real, 2}((2,3))
        v = Vector{Real}(undef, arrshape)
        @test @inferred(length(v)) == 6
        @test @inferred(size(v)) == (6,)

        data1 = [1;2;3;4;7;8;9]
        scalarshape = ScalarShape{Real}()
        ntshape = NamedTupleShape(a=arrshape, b=scalarshape)
        shapedasnt = ntshape(data1)
        @test stripscalar(Ref(shapedasnt)) == Ref(shapedasnt)[]

        @test_throws ArgumentError Broadcast.broadcastable(ntshape)

        named_shapes = (
            a = ArrayShape{Real}(2, 3),
            b = ScalarShape{Real}(),
            c = ConstValueShape(4.2),
            x = ConstValueShape([11 21; 12 22]),
            y = ArrayShape{Real}(4)
        )
        shape = NamedTupleShape(;named_shapes...)
        sntshape = NamedTupleShape(ShapedAsNT; named_shapes...)
        data2 = VectorOfSimilarVectors(reshape(collect(1:22), 11, 2))
        @test_throws ArgumentError ValueShapes._checkcompat_inner(ntshape, data2)
        @test ValueShapes._checkcompat_inner(shape, data2) == nothing

        let vs = sntshape, x = rand(totalndof(vs)), xs = nestedview(rand(totalndof(vs), 5))
            vs_jacobian(f, x) = 1

            InverseFunctions.test_inverse(vs, x)
            ChangesOfVariables.test_with_logabsdet_jacobian(vs, x, vs_jacobian)
            ChangesOfVariables.test_with_logabsdet_jacobian(inverse(vs), vs(x), vs_jacobian)

            bc_vs = Base.Fix1(broadcast, vs)
            bc_unshaped = Base.Fix1(broadcast, unshaped)
            InverseFunctions.test_inverse(bc_vs, xs)
            @test with_logabsdet_jacobian(bc_vs, xs)[1] isa ShapedAsNTArray
            ChangesOfVariables.test_with_logabsdet_jacobian(bc_vs, xs, vs_jacobian)
            @test with_logabsdet_jacobian(inverse(bc_vs), vs.(xs))[1] isa ArrayOfSimilarArrays
            ChangesOfVariables.test_with_logabsdet_jacobian(inverse(bc_vs), vs.(xs), vs_jacobian)
            @test with_logabsdet_jacobian(bc_unshaped, vs.(xs))[1] isa ArrayOfSimilarArrays
            ChangesOfVariables.test_with_logabsdet_jacobian(bc_unshaped, vs.(xs), vs_jacobian)

            for f in [vs, inverse(vs), bc_vs, inverse(bc_vs), unshaped]
                @test @inferred(identity ∘ f) === f
                @test @inferred(f ∘ identity) === f
            end
        end

        @test @inferred(unshaped(4.2)) isa Fill{Float64,1}
        @test unshaped(4.2) == [4.2]
        @test @inferred(unshaped(view([4.2], 1))) isa SubArray{Float64,1,Vector{Float64}}
        @test unshaped(view([4.2], 1)) == [4.2]
        @test @inferred(unshaped(Array(view([4.2], 1)))) isa SubArray{Float64,1,Vector{Float64}}
        @test unshaped(Array(view([4.2], 1))) == [4.2]
        let x = rand(15)
            @test @inferred(unshaped(x)) === x
            @test @inferred(unshaped(Base.ReshapedArray(x, (3, 5), ()))) === x

            @test @inferred(broadcast(unshaped, x)) isa ArrayOfSimilarArrays{Float64,1,1,2,<:Base.ReshapedArray}
            @test broadcast(unshaped, x) == nestedview(reshape(x, 1, 15))
        end

        let A = rand(1,15)
            @test @inferred(broadcast(unshaped, view(A, 1, :))) isa ArrayOfSimilarArrays{Float64,1,1,2,<:SubArray}
            @test broadcast(unshaped, view(A, 1, :)) == nestedview(A)
        end
    end

    @testset "realnumtype" begin
        a = 4.2f0
        b = Complex(a, a)
        A = fill(fill(rand(Float32, 5), 10), 5)
        B = fill(rand(Float32, 5), 10)
        C = ArrayOfSimilarArrays(B)
        nt = (a = 4.2f0, b = 42)
        tpl = (4.2f0, Float16(2.3))
        deepnt = (a = a, b = b, A = A, B = B, C = C, nt = nt, tpl = tpl)
        for x in [a, A, B, C, nt, tpl, deepnt]
            @test @inferred (realnumtype(typeof(x))) == Float32
        end
    end

    @testset "value shape comparison" begin
        for (T, U) in [(Real, Float64), (Float64, Real), (Real, Real), (Float64, Float64)]
            @test @inferred(ScalarShape{T}() <= ScalarShape{U}())[1] == (T <: U)
            @test @inferred(rrule(Base.:(<=), ScalarShape{T}(), ScalarShape{U}()))[1] == (T <: U)
            @test @inferred(rrule(Base.:(<=), ScalarShape{T}(), ScalarShape{U}()))[2](T <: U) == (NoTangent(), NoTangent(), NoTangent())

            @test @inferred(ScalarShape{T}() >= ScalarShape{U}())[1] == (T >: U)
            @test @inferred(rrule(Base.:(>=), ScalarShape{T}(), ScalarShape{U}()))[1] == (T >: U)
            @test @inferred(rrule(Base.:(>=), ScalarShape{T}(), ScalarShape{U}()))[2](T >: U) == (NoTangent(), NoTangent(), NoTangent())

            @test @inferred(ScalarShape{T}() == ScalarShape{U}())[1] == (T == U)
            @test @inferred(rrule(Base.:(==), ScalarShape{T}(), ScalarShape{U}()))[1] == (T == U)
            @test @inferred(rrule(Base.:(==), ScalarShape{T}(), ScalarShape{U}()))[2](T == U) == (NoTangent(), NoTangent(), NoTangent())
        end
    end
end


text:
using HNSW
import HNSW: LayeredGraph, add_vertex!, add_edge!, get_entry_level, levelof, neighbors,max_connections
using Test
using LinearAlgebra
using NearestNeighbors

@testset "Nearest & Furthest" begin
    for i=1:10
        data = [rand(10) for i = 1:100]
        q = rand(10)
        nearest_idx = findmin([norm(q-x) for x ∈ data])[2]
        furthest_idx = findmax([norm(q-x) for x ∈ data])[2]

        #prepare candidates
        cand = HNSW.NeighborSet(1, norm(q-data[1]))
        for n = 2:100
            HNSW.insert!(cand, HNSW.Neighbor(n, norm(q-data[n])))
        end
        #nearest
        @test nearest_idx == HNSW.nearest(cand).idx

        #extract_nearest!
        @test nearest_idx == HNSW.pop_nearest!(cand).idx
        @test nearest_idx != HNSW.nearest(cand)

        #furthest
        @test furthest_idx == HNSW.furthest(cand).idx
        HNSW.pop_furthest!(cand)
        @test furthest_idx != HNSW.furthest(cand)
    end
end

#Test Visited List Pool
@testset "VisitedListPool" begin
    elnum = 20
    vlp = HNSW.VisitedListPool(1,elnum)
    vl = HNSW.get_list(vlp)
    println(vl)
    for i = 1:elnum
        @test HNSW.isvisited(vl, i) == false
        HNSW.visit!(vl,i)
        @test HNSW.isvisited(vl, i) == true
    end
    for i=1:500
         HNSW.reset!(vl)
         j = rand(1:elnum)
         @test HNSW.isvisited(vl, j) == false
    end
    println(vl)
end

#Test layered Graph
@testset "LayeredGraph" begin
    num_elements = 100
    M, M0 = 16, 32
    m_L = 1.
    lg = HNSW.LayeredGraph{UInt32}(num_elements, M0, M, m_L)
    @test max_connections(lg, 1) == M0
    @test max_connections(lg, 2) == M
    @testset "add_vertex! & get_entry_level" for i=1:10
        level = rand(1:10)
        add_vertex!(lg, i, level)
        #@test get_entry(lg) >= level
    end
    @testset "add_edge!" begin
        for i = 1:10, j=1:10
            if i==j
                #@test_throws AssertionError add_edge!(lg, 1, i, j)
            else
                level = rand(1:5)
                if levelof(lg, i) < level || levelof(lg, j) < level
                    #@test_throws AssertionError add_edge!(lg, level, i,j)
                else
                    add_edge!(lg, level, i,j)

                    @test j ∈ [neighbors(lg,level,i)...]
                end
            end
        end
    end
    @testset "neighbors" begin
        for i = 1:10
            level = rand(1:10)
            if level > levelof(lg, i)
                #@test_throws AssertionError neighbors(lg, i, level)
            else
                #retrieve neighbors
                idx_offset = level > 1 ? lg.M0 + lg.M*(level-2) : 0
                maxM = HNSW.max_connections(lg, level)
                links = lg.linklist[i][idx_offset.+(1:maxM)]
                links = links[links .> 0]
                @test [neighbors(lg, level, i)...] == links
            end
        end
    end
end


text:
function addto!(jm::ModelView, m::EnergyModel, sn::SubNetwork)
    T = axis(sn.model, :snapshots)
    branches = sn[PassiveBranch]
    B = axis(branches)
    length(B) == 0 && return

    p = get(branches, :p, B, T)
    effimp = mapcat(effectiveimpedance, branches)

    C = cycle_matrix(sn)
    Cl = rowvals(C)
    Cv = nonzeros(C)

    @constraint(jm, cycles[c=1:size(C,2),t=T], sum(Cv[j] * effimp[Cl[j]] * p[Cl[j],t] for j=nzrange(C,c)) == 0)
end

addcomponent(SubNetwork, :subnetworks)

max_stars_repo_name: coroa/EnergyModels

max_issues_repo_name: coroa/EnergyModels

max_forks_repo_name: coroa/EnergyModels


text:
# Activate project environment
import Pkg; Pkg.activate(".")

Pkg.instantiate()

# Load required packages
using CSV
using DataFrames
using Distributions
using EcologicalNetworks
import Mangal
using Plots
using Random
using StatsPlots
using StatsBase

# Load functions
include("src/predict_L.jl")
include("src/maxent_dd.jl")
include("src/empirical_dd.jl")

# Make figures
include("scripts/00_import_mangal_data.jl")
include("scripts/01_flexible_links_posterior.jl")
include("scripts/02_dd_speculative_webs.jl")
include("scripts/03_dd_mangal_webs.jl")

max_stars_repo_name: FrancisBanville/Trophic-METE

max_issues_repo_name: FrancisBanville/Trophic-METE

max_forks_repo_name: FrancisBanville/Trophic-METE


text:
using FIRITools
using Test, Statistics, DelimitedFiles
using TypedTables

@testset "FIRITools.jl" begin
    ALTITUDE_LENGTH = 89
    @test size(FIRITools.DATA) == (ALTITUDE_LENGTH, 1980)  # for firi2018
    @test length(FIRITools.HEADER) == 1980
    @test length(columnnames(FIRITools.HEADER)) == 5

    @test FIRITools.values(:f10_7) == FIRITools.values("f10_7")
    @test FIRITools.values(:f10_7) == unique(FIRITools.HEADER.f10_7)
    @test FIRITools.values(:chi) == unique(FIRITools.HEADER.chi)
    @test FIRITools.values(:lat) == unique(FIRITools.HEADER.lat)
    @test FIRITools.values(:month) == unique(FIRITools.HEADER.month)

    @test FIRITools.twoclosest(FIRITools.HEADER.chi, 90) == [90, 90]
    @test FIRITools.twoclosest(FIRITools.HEADER.chi, 89) == [85, 90]
    @test FIRITools.twoclosest(FIRITools.HEADER.chi, 135) == [100, 130]

    @test size(FIRITools.selectprofiles()) == (ALTITUDE_LENGTH, 1980)
    @test size(FIRITools.selectprofiles(month=(1, 6))) == (ALTITUDE_LENGTH, 1980÷2)
    @test_logs (:warn, "15 is not in model parameters. Looking for interpolating form of `selectprofiles`?") FIRITools.selectprofiles(chi=15)

    # Try interpolating chi and lat. Should always have same size
    @test size(FIRITools.selectprofiles(15, 10)) == (ALTITUDE_LENGTH, 36)  # interpolate both
    @test size(FIRITools.selectprofiles(15.1, 82.3)) == (ALTITUDE_LENGTH, 36)  # extrapolate
    @test size(FIRITools.selectprofiles(30, 30)) == (ALTITUDE_LENGTH, 36)  # no interpolation
    @test size(FIRITools.selectprofiles(30, 10)) == (ALTITUDE_LENGTH, 36) # interpolate lat
    @test size(FIRITools.selectprofiles(29, 30)) == (ALTITUDE_LENGTH, 36) # interpolate chi

    profs85 = FIRITools.selectprofiles(85, 31)
    profs89 = FIRITools.selectprofiles(89, 31)
    profs90 = FIRITools.selectprofiles(90, 31)
    @test all(profs85 .> profs89 .> profs90)
    # This test doesn't work for lat because the profiles are largely overlapping.
    # Even using `firi` instead of `selectprofiles` doesn't work because the mean of the
    # profiles causes there to be some crossing density regions

    @test size(firi()) == (ALTITUDE_LENGTH,)
    @test size(firi(month=(1, 6))) == (ALTITUDE_LENGTH,)
    @test size(firi(15, 10)) == (ALTITUDE_LENGTH,)

    @test_throws ArgumentError firi(-10, 45)
    @test_throws ArgumentError firi(30, -15)
    @test_logs (:warn, "`chi` greater than 130° uses `chi = 130°`") firi(160, 45)

    @test FIRITools.quantile(1, chi=(0, 90)) ==
        Statistics.quantile.(eachrow(FIRITools.selectprofiles(chi=(0, 90))), 1)

    q1 = FIRITools.quantile(chi=(0, 90), (0.3, 0.5))
    @test q1[:,1] == FIRITools.quantile(chi=(0, 90), 0.3)
    @test q1[:,2] == FIRITools.quantile(chi=(0, 90), 0.5)

    function quantilewoextrap(ref, p)
        # Compare to reference profiles without extrapolation
        v, h = ref  # values, header
        refmask = [a in FIRITools.ALTITUDE[:,1] for a in v[:,1]]
        fmask = [a in v[:,1] for a in FIRITools.ALTITUDE[:,1]]
        isapprox(v[refmask,2],
                 FIRITools.quantile(chi=(0, 95), p)[fmask], rtol=0.01)
    end

    @test quantilewoextrap(readdlm("FIRIday_30.csv", ','; header=true), 0.3)
    @test quantilewoextrap(readdlm("FIRIday_90.csv", ','; header=true), 0.9)

    # plot(ref30.ne, ref30.alt, label="Wei 30", xlims=(10^4, 4e10), ylims=(55, 110), xscale=:log10)
    # plot!(FIRITools.quantile(chi=(0, 95), 0.3), FIRITools.ALTITUDE, label=30, xscale=:log10)

    firialt = minimum(FIRITools.ALTITUDE):1000:maximum(FIRITools.ALTITUDE)
    p1 = FIRITools.extrapolate(firialt, firi(), 30e3:1e3:120e3)
    p3 = FIRITools.extrapolate(firi(), 30e3:1e3:120e3)
    p5 = FIRITools.extrapolate(firi(), 70e3:1e3:120e3)

    @test !any(isinf, p1)
    @test !any(isnan, p1)

    @test length(p1) == length(30:120)
    @test length(p5) == length(70:120)

    @test p1 == p3
    @test p1[65 .< 30:120] ≈ firi()[65e3 .< FIRITools.ALTITUDE .<= 120e3]
    @test p1[70 .<= 30:120] == p5

    # plot(firi(), FIRITools.ALTITUDE/1000, xscale=:log10, legend=false)
    # plot!(FIRITools.extrapolate(firi(), (30:148)*1e3; N=6), 30:148)
    # plot!(p5, 70:120)

    # ep = FIRITools.extrapolate(FIRITools.quantile(chi=(0, 95), 0.5), 40:110; max_altitude=60)
    # plot(ep, 40:110, xscale=:log10)

    p11 = FIRITools.extrapolate(FIRITools.DATA[:,[1, 2]], 30e3:1e3:120e3)
    @test p11[:,1] == FIRITools.extrapolate(FIRITools.DATA[:,1], 30e3:1e3:120e3)
    @test p11[:,2] == FIRITools.extrapolate(FIRITools.DATA[:,2], 30e3:1e3:120e3)

    newz = 30e3:250:110e3
    @test length(FIRITools.extrapolate(firi(), newz)) == length(newz)
end


text:
struct polytope
    A;
    b;
    l;
    u;
    model;
    x;
end

function LP(c,P)

    @objective(P.model, Min, c'* P.x)
    optimize!(P.model)
    return value.(P.x)

end

function randomDoublyStochasticMatrix(n; num_perm = n^2)

    M = zeros(n,n);
    α = rand(num_perm)

    α = α / sum(α);

    for i=1:num_perm
        perm = randperm(n);
        for j=1:n
            M[perm[j],j] += α[i];
        end
    end

    return M;
end

function birkhoffPolytope(n)

    M = zeros(n*n,2*n);
    for i=1:n
        M[(i-1)*n*n + (i-1)*n + 1 : (i-1)*n*n + (i-1)*n + n ] = ones(n,1);
    end
    for i=1:n
        for j=1:n
            M[n*n*n + (i-1)*n*n + (j-1)*n + i] = 1;
        end
    end

    A = sparse(M');
    b = ones(2*n);

    model = Model(with_optimizer(Clp.Optimizer))
    set_silent(model)
    @variable(model, 0 <= x[1:n*n] <= 1)
    @constraint(model, A * x .== b)

    return polytope(A,b,0,1,model,x)

end


text:
#!/usr/bin/julia

using DataFrames

srand(46)

SAMPLES = 100
ITERS=200
ITERS_p=200

values = readdlm("SDvalues.dat") |> vec
# A = 5 + randn(10000)
# B = -5 + randn(10000)
# values = vcat(A,B)
info("Real mean: $(mean(values))")

function bootstrap(array, frac; N=1000)
    L = length(array)
    chunksize = floor(Int64,frac*L)
    0 < chunksize ≤ L || error("chunksize = $chunksize ∉ (0,L)")
    means = zeros(N)
    for i ∈ 1:N
        chunk = rand(array,chunksize)
        means[i] = mean(chunk)
    end
    return mean(means), std(means)
end

function jacknife(array)
    L = length(array)
    means = zeros(L)
    for i ∈ 1:L
        means[i] = mean(array[j] for j in 1:L if j≠i)
    end
    return mean(means), std(means)
end

jdf = DataFrame()
bdf = DataFrame()

fr = 1/e
for i in 1:ITERS
    smallsample = rand(values,SAMPLES)
    jackμ, jackσ = jacknife(smallsample)
    bootμ, bootσ = bootstrap(smallsample, fr)
    jdf = vcat(jdf, DataFrame(i=i,
                              μ=jackμ,
                              σ=jackσ))
    bdf = vcat(bdf, DataFrame(i=i,
                              fraction=fr,
                              μ=bootμ,
                              σ=bootσ))
end

writetable("data_jack.csv", jdf)
writetable("data_boot.csv", bdf)

max_stars_repo_name: RedPointyJackson/tfg

max_issues_repo_name: RedPointyJackson/tfg

max_forks_repo_name: RedPointyJackson/tfg


text:
@testset "val" begin
    x = RealSemiringElement(3.0)
    @test val(x) ≈ 3
    @test typeof(val(x)) == eltype(x)
end


text:
module ExitWithParent

export exit_with_parent!

"""
    exit_with_parent!(p::Base.Process, signal = 2)

When this julia process exits, the provided process `p` will be killed with `signal`.

This spawns a child process that simply polls until the parent is killed, and then kills `p`.
"""
function exit_with_parent! end
if Sys.isunix()
    function exit_with_parent!(p::Base.Process, signal = 2)  # SIGINT (ctrl-c)
        targetpid = getpid(p)
        orig_parent = getpid()
        open(`$(Base.julia_cmd()) -e "
            while true
                ppid = ccall(:getppid, Cint, ())
                sleep(1)
                if ppid != $orig_parent
                    ccall(:kill, Cint, (Cint,Cint), $targetpid, $signal)
                    exit()
                end
            end
        "`);
        nothing
    end
else
    # TODO: What should this do for other systems?
    function exit_with_parent!(p::Base.Process, signal = 2) end # SIGINT (ctrl-c)
end

end # module

max_stars_repo_name: NHDaly/ExitWithParent

max_issues_repo_name: NHDaly/ExitWithParent

max_forks_repo_name: NHDaly/ExitWithParent


text:
using FluxGAN
using MLDatasets
using Flux
using CairoMakie
using LaTeXStrings

const n = parse(Int, ARGS[1])
const m = parse(Int, ARGS[2])
const skip = parse(Int, ARGS[3])

images = MNIST.traintensor(Float32)

const img_size = size(images)[1:end-1]
const img_num  = size(images)[end]
const img_dim  = *(img_size...)

# centering pixels, mapping to [-1,1]; and flattening image tensors
images = reshape((@. 2f0 * images - 1f0), (img_dim, img_num))

const d = 100

generator = Chain(
    Dense(d, 1200, relu),
    Dense(1200, 1200, relu),
    Dense(1200, img_dim, tanh_fast)
)

discriminator = Chain(
    Maxout(() -> Dense(img_dim, 240), 5),
    Maxout(() -> Dense(240, 240), 5),
    Dense(240, 1)
)

# minibatch size

ks = [1, 5, 20, 40]
titles = [latexstring("\$k = \$$k") for k in ks ]

fig = Figure()
ax1 = Axis(fig[1,1], xlabel="iteration", ylabel="loss", title=titles[1])
ax2 = Axis(fig[1,2], xlabel="iteration", ylabel="loss", title=titles[2])
ax3 = Axis(fig[2,1], xlabel="iteration", ylabel="loss", title=titles[3])
ax4 = Axis(fig[2,2], xlabel="iteration", ylabel="loss", title=titles[4])
axs = [ax1, ax2, ax3, ax4]

ts = collect(skip:skip:n)

maxs = []

for (i, k) in enumerate(ks)
    model = GAN(deepcopy(generator), deepcopy(discriminator);
        minibatch=m,
        latent_dim=d,
        img_size=img_size,
        dscr_loops=k
    )
    gen_losses, dscr_losses = train!(model, images;
                                     iterations=n,
                                     skip=skip,
                                     return_losses=true)
    push!(maxs, maximum([gen_losses; dscr_losses]))
    lines!(axs[i], ts, gen_losses; label=L"G", color=:green)
    lines!(axs[i], ts, dscr_losses; label=L"D", color=:purple)
    axislegend(axs[i])
end

max = maximum(maxs)

for ax in axs
    ylims!(ax, 0, max)
end

save("images/loss_plots/loss_for_ks_n_$n.png", fig)

max_stars_repo_name: aarontrowbridge/GAN

max_issues_repo_name: aarontrowbridge/GAN

max_forks_repo_name: aarontrowbridge/GAN


text:
using DelimitedFiles
function x(mv)
	mp = readdlm("3.txt")
	pos = [1,1]
	n = 0
	while pos[2] <= length(mp)
		if mp[pos[2]][pos[1]] == '#'
			n += 1
		end
		pos += mv
		pos[1] = (pos[1]-1) % length(mp[1]) + 1
	end
	return n
end
# part 1
println(x([3,1]))
# part 2
println(x([1,1])*x([3,1])*x([5,1])*x([7,1])*x([1,2]))


text:
@testset "Cubic equation" begin
    Base.info("Testing cubic equation...")
    solver = Solver()
    @test isapprox(solver(1.0,-3.0,3.0,-1.0), ones(3))
    @test isapprox(solver(2.0,-3.0,-5.0,6.0), [-1.5,1.0,2.0])
    @test isapprox.(solver(1.0,1.0,1.0,-3.0), Union{Real,Complex}[1.0,-1.0-im*√(2),-1.0+im*√(2)]) == ones(Bool, 3)
    @test_throws ErrorException solver(0.0,1.0,1.0,-3.0)
    Base.info("Done testing cubic equation.")
end


text:
module DataStructures

    import Base: <, <=, ==, length, isempty, iterate,
                 show, dump, empty!, getindex, setindex!, get, get!,
                 in, haskey, keys, merge, copy, cat, collect,
                 push!, pop!, pushfirst!, popfirst!, insert!, lastindex,
                 union!, delete!, similar, sizehint!, empty, append!,
                 isequal, hash, map, filter, reverse, peek,
                 first, last, eltype, getkey, values, sum,
                 merge, merge!, lt, Ordering, ForwardOrdering, Forward,
                 ReverseOrdering, Reverse, Lt,
                 isless, union, intersect, symdiff, setdiff, issubset,
                 searchsortedfirst, searchsortedlast, in,
                 eachindex, keytype, valtype, minimum, maximum, size,
                 zero, checkbounds

    using OrderedCollections
    import OrderedCollections: filter, filter!, isordered
    export OrderedDict, OrderedSet, LittleDict
    export DefaultDict, DefaultOrderedDict

    export complement, complement!

    export Deque, Stack, Queue, CircularDeque
    export enqueue!, dequeue!, dequeue_pair!, update!, reverse_iter
    export capacity, num_blocks, top_with_handle, sizehint!

    export Accumulator, counter, reset!, inc!, dec!
    export IntDisjointSets, DisjointSets, num_groups, find_root!, in_same_set, root_union!
    export FenwickTree, length, inc!, dec!, incdec!, prefixsum

    export AbstractHeap, compare, extract_all!
    export BinaryHeap, BinaryMinHeap, BinaryMaxHeap, nlargest, nsmallest
    export MutableBinaryHeap, MutableBinaryMinHeap, MutableBinaryMaxHeap
    export heapify!, heapify, heappop!, heappush!, isheap
    export BinaryMinMaxHeap, popmin!, popmax!, popall!

    export Trie, subtrie, keys_with_prefix, partial_path

    export LinkedList, Nil, Cons, nil, cons, head, tail, list, filter, cat,
           reverse
    export MutableLinkedList
    export SortedDict, SortedMultiDict, SortedSet
    export SDToken, SDSemiToken, SMDToken, SMDSemiToken
    export SetToken, SetSemiToken
    export startof
    export pastendsemitoken, beforestartsemitoken
    export searchsortedafter, searchequalrange
    export packcopy, packdeepcopy
    export exclusive, inclusive, semitokens
    export orderobject, ordtype, Lt, compare, onlysemitokens

    export MultiDict, enumerateall
    export RobinDict

    export DiBitVector

    export findkey

    include("delegate.jl")

    include("deque.jl")
    include("circ_deque.jl")
    include("stack.jl")
    include("queue.jl")
    include("accumulator.jl")
    include("disjoint_set.jl")
    include("heaps.jl")

    include("default_dict.jl")
    include("dict_support.jl")
    include("trie.jl")

    include("int_set.jl")

    include("fenwick.jl")

    include("list.jl")
    include("mutable_list.jl")
    include("balanced_tree.jl")
    include("tokens.jl")

    import .Tokens: IntSemiToken

    include("multi_dict.jl")
    include("sorted_dict.jl")
    include("sorted_multi_dict.jl")
    include("sorted_set.jl")
    include("tokens2.jl")
    include("container_loops.jl")
    include("robin_dict.jl")
    export
        CircularBuffer,
        capacity,
        isfull
    include("circular_buffer.jl")

    export status
    export deref_key, deref_value, deref, advance, regress

    export PriorityQueue, peek

    include("priorityqueue.jl")
    include("sparse_int_set.jl")
    export SparseIntSet

    include("dibit_vector.jl")
    include("deprecations.jl")
end