text:
module GLU
#Autogenerated, hopefully stays that way!

import GetC.@getCFun

const libGLU = @windows? "glu32" : "libGLU"

typealias GLenum Cuint
typealias GLboolean Cuchar
typealias GLbitfield Cuint
typealias GLvoid Void
typealias GLbyte Cuchar
typealias GLshort Cshort
typealias GLint Cint
typealias GLubyte Cuchar
typealias GLushort Cushort
typealias GLuint Cuint
typealias GLsizei Cint
typealias GLfloat Cfloat
typealias GLclampf Cfloat
typealias GLdouble Cdouble
typealias GLclampd Cdouble
typealias ptrdiff_t Clonglong
typealias size_t Clonglong
typealias wchar_t Cint
typealias GLchar Cuchar
typealias GLintptr ptrdiff_t
typealias GLsizeiptr ptrdiff_t
typealias GLintptrARB ptrdiff_t
typealias GLsizeiptrARB ptrdiff_t
typealias GLcharARB Cuchar
typealias GLhandleARB Cuint
typealias GLhalfARB Cushort
typealias GLhalfNV Cushort
typealias int8_t Cuchar
typealias int16_t Cshort
typealias int32_t Cint
typealias int64_t Clonglong
typealias uint8_t Cuchar
typealias uint16_t Cushort
typealias uint32_t Cuint
typealias uint64_t Culonglong
typealias int_least8_t Cuchar
typealias int_least16_t Cshort
typealias int_least32_t Cint
typealias int_least64_t Clonglong
typealias uint_least8_t Cuchar
typealias uint_least16_t Cushort
typealias uint_least32_t Cuint
typealias uint_least64_t Culonglong
typealias int_fast8_t Cuchar
typealias int_fast16_t Clonglong
typealias int_fast32_t Clonglong
typealias int_fast64_t Clonglong
typealias uint_fast8_t Cuchar
typealias uint_fast16_t Culonglong
typealias uint_fast32_t Culonglong
typealias uint_fast64_t Culonglong
typealias intptr_t Clonglong
typealias uintptr_t Culonglong
typealias intmax_t Clonglong
typealias uintmax_t Culonglong
typealias __gwchar_t Cint
export GLenum
export GLboolean
export GLbitfield
export GLvoid
export GLbyte
export GLshort
export GLint
export GLubyte
export GLushort
export GLuint
export GLsizei
export GLfloat
export GLclampf
export GLdouble
export GLclampd
export ptrdiff_t
export size_t
export wchar_t
export GLchar
export GLintptr
export GLsizeiptr
export GLintptrARB
export GLsizeiptrARB
export GLcharARB
export GLhandleARBt32
export GLhalfARB
export GLhalfNV
export int8_t
export int16_t
export int32_t
export int64_t
export uint8_t
export uint16_t
export uint32_t
export uint64_t
export int_least8_t
export int_least16_t
export int_least32_t
export int_least64_t
export uint_least8_t
export uint_least16_t
export uint_least32_t
export uint_least64_t
export int_fast8_t
export int_fast16_t
export int_fast32_t
export int_fast64_t
export uint_fast8_t
export uint_fast16_t
export uint_fast32_t
export uint_fast64_t
export intptr_t
export uintptr_t
export intmax_t
export uintmax_t
export __gwchar_t
@getCFun libGLU imaxabs imaxabs(__n::intmax_t)::intmax_t
export imaxabs
@getCFun libGLU imaxdiv imaxdiv(__numer::intmax_t,__denom::intmax_t)::imaxdiv_t
export imaxdiv
@getCFun libGLU strtoimax strtoimax(__nptr::Ptr,__endptr::Ptr,__base::Int32)::intmax_t
export strtoimax
@getCFun libGLU strtoumax strtoumax(__nptr::Ptr,__endptr::Ptr,__base::Int32)::uintmax_t
export strtoumax
@getCFun libGLU wcstoimax wcstoimax(__nptr::Ptr,__endptr::Ptr,__base::Int32)::intmax_t
export wcstoimax
@getCFun libGLU wcstoumax wcstoumax(__nptr::Ptr,__endptr::Ptr,__base::Int32)::uintmax_t
export wcstoumax
typealias GLint64EXT int64_t
typealias GLuint64EXT uint64_t
typealias GLint64 int64_t
typealias GLuint64 uint64_t
typealias GLvdpauSurfaceNV GLintptr
typealias GLeglImageOES Ptr
export GLint64EXT
export GLuint64EXT
export GLint64
export GLuint64
export GLvdpauSurfaceNV
export GLeglImageOES
@getCFun libGLU gluPerspective gluPerspective(fovy::GLdouble,aspect::GLdouble,znear::GLdouble,zfar::GLdouble)::Void
export gluPerspective
@getCFun libGLU gluBuild2DMipmaps gluBuild2DMipmaps(target::GLenum,internalFormat::GLint,width::GLsizei,height::GLsizei,format::GLenum,thetype::GLenum,data::Ptr{Uint8})::GLint
export gluBuild2DMipmaps
@getCFun libGLU gluNewQuadric gluNewQuadric()::Ptr{Void}
export gluNewQuadric
@getCFun libGLU gluQuadricNormals gluQuadricNormals(quad::Ptr{Void},normal::GLenum)::Void
export gluQuadricNormals
@getCFun libGLU gluQuadricTexture gluQuadricTexture(quad::Ptr{Void},texture::GLboolean)::Void
export gluQuadricTexture
@getCFun libGLU gluCylinder gluCylinder(quad::Ptr{Void},base::GLdouble,top::GLdouble,height::GLdouble,slices::GLint,stacks::GLint)::Void
export gluCylinder
@getCFun libGLU gluDisk gluDisk(quad::Ptr{Void},inner::GLdouble,outer::GLdouble,slices::GLint,stacks::GLint)::Void
export gluDisk
@getCFun libGLU gluSphere gluSphere(quad::Ptr{Void},radius::GLdouble,slices::GLint,stacks::GLint)::Void
export gluSphere
@getCFun libGLU gluPartialDisk gluPartialDisk(quad::Ptr{Void},inner::GLdouble,outer::GLdouble,slices::GLint,loops::GLint,start::GLdouble,sweep::GLdouble)::Void
export gluPartialDisk
const GLU_EXT_object_space_tess       = 1
const GLU_EXT_nurbs_tessellator       = 1
const GLU_FALSE                       = 0
const GLU_TRUE                        = 1
const GLU_VERSION_1_1                 = 1
const GLU_VERSION_1_2                 = 1
const GLU_VERSION_1_3                 = 1
const GLU_VERSION                     = 100800
const GLU_EXTENSIONS                  = 100801
const GLU_INVALID_ENUM                = 100900
const GLU_INVALID_VALUE               = 100901
const GLU_OUT_OF_MEMORY               = 100902
const GLU_INCOMPATIBLE_GL_VERSION     = 100903
const GLU_INVALID_OPERATION           = 100904
const GLU_OUTLINE_POLYGON             = 100240
const GLU_OUTLINE_PATCH               = 100241
const GLU_NURBS_ERROR                 = 100103
const GLU_ERROR                       = 100103
const GLU_NURBS_BEGIN                 = 100164
const GLU_NURBS_BEGIN_EXT             = 100164
const GLU_NURBS_VERTEX                = 100165
const GLU_NURBS_VERTEX_EXT            = 100165
const GLU_NURBS_NORMAL                = 100166
const GLU_NURBS_NORMAL_EXT            = 100166
const GLU_NURBS_COLOR                 = 100167
const GLU_NURBS_COLOR_EXT             = 100167
const GLU_NURBS_TEXTURE_COORD         = 100168
const GLU_NURBS_TEX_COORD_EXT         = 100168
const GLU_NURBS_END                   = 100169
const GLU_NURBS_END_EXT               = 100169
const GLU_NURBS_BEGIN_DATA            = 100170
const GLU_NURBS_BEGIN_DATA_EXT        = 100170
const GLU_NURBS_VERTEX_DATA           = 100171
const GLU_NURBS_VERTEX_DATA_EXT       = 100171
const GLU_NURBS_NORMAL_DATA           = 100172
const GLU_NURBS_NORMAL_DATA_EXT       = 100172
const GLU_NURBS_COLOR_DATA            = 100173
const GLU_NURBS_COLOR_DATA_EXT        = 100173
const GLU_NURBS_TEXTURE_COORD_DATA    = 100174
const GLU_NURBS_TEX_COORD_DATA_EXT    = 100174
const GLU_NURBS_END_DATA              = 100175
const GLU_NURBS_END_DATA_EXT          = 100175
const GLU_NURBS_ERROR1                = 100251
const GLU_NURBS_ERROR2                = 100252
const GLU_NURBS_ERROR3                = 100253
const GLU_NURBS_ERROR4                = 100254
const GLU_NURBS_ERROR5                = 100255
const GLU_NURBS_ERROR6                = 100256
const GLU_NURBS_ERROR7                = 100257
const GLU_NURBS_ERROR8                = 100258
const GLU_NURBS_ERROR9                = 100259
const GLU_NURBS_ERROR10               = 100260
const GLU_NURBS_ERROR11               = 100261
const GLU_NURBS_ERROR12               = 100262
const GLU_NURBS_ERROR13               = 100263
const GLU_NURBS_ERROR14               = 100264
const GLU_NURBS_ERROR15               = 100265
const GLU_NURBS_ERROR16               = 100266
const GLU_NURBS_ERROR17               = 100267
const GLU_NURBS_ERROR18               = 100268
const GLU_NURBS_ERROR19               = 100269
const GLU_NURBS_ERROR20               = 100270
const GLU_NURBS_ERROR21               = 100271
const GLU_NURBS_ERROR22               = 100272
const GLU_NURBS_ERROR23               = 100273
const GLU_NURBS_ERROR24               = 100274
const GLU_NURBS_ERROR25               = 100275
const GLU_NURBS_ERROR26               = 100276
const GLU_NURBS_ERROR27               = 100277
const GLU_NURBS_ERROR28               = 100278
const GLU_NURBS_ERROR29               = 100279
const GLU_NURBS_ERROR30               = 100280
const GLU_NURBS_ERROR31               = 100281
const GLU_NURBS_ERROR32               = 100282
const GLU_NURBS_ERROR33               = 100283
const GLU_NURBS_ERROR34               = 100284
const GLU_NURBS_ERROR35               = 100285
const GLU_NURBS_ERROR36               = 100286
const GLU_NURBS_ERROR37               = 100287
const GLU_AUTO_LOAD_MATRIX            = 100200
const GLU_CULLING                     = 100201
const GLU_SAMPLING_TOLERANCE          = 100203
const GLU_DISPLAY_MODE                = 100204
const GLU_PARAMETRIC_TOLERANCE        = 100202
const GLU_SAMPLING_METHOD             = 100205
const GLU_U_STEP                      = 100206
const GLU_V_STEP                      = 100207
const GLU_NURBS_MODE                  = 100160
const GLU_NURBS_MODE_EXT              = 100160
const GLU_NURBS_TESSELLATOR           = 100161
const GLU_NURBS_TESSELLATOR_EXT       = 100161
const GLU_NURBS_RENDERER              = 100162
const GLU_NURBS_RENDERER_EXT          = 100162
const GLU_OBJECT_PARAMETRIC_ERROR     = 100208
const GLU_OBJECT_PARAMETRIC_ERROR_EXT = 100208
const GLU_OBJECT_PATH_LENGTH          = 100209
const GLU_OBJECT_PATH_LENGTH_EXT      = 100209
const GLU_PATH_LENGTH                 = 100215
const GLU_PARAMETRIC_ERROR            = 100216
const GLU_DOMAIN_DISTANCE             = 100217
const GLU_MAP1_TRIM_2                 = 100210
const GLU_MAP1_TRIM_3                 = 100211
const GLU_POINT                       = 100010
const GLU_LINE                        = 100011
const GLU_FILL                        = 100012
const GLU_SILHOUETTE                  = 100013
const GLU_SMOOTH                      = 100000
const GLU_FLAT                        = 100001
const GLU_NONE                        = 100002
const GLU_OUTSIDE                     = 100020
const GLU_INSIDE                      = 100021
const GLU_TESS_BEGIN                  = 100100
const GLU_BEGIN                       = 100100
const GLU_TESS_VERTEX                 = 100101
const GLU_VERTEX                      = 100101
const GLU_TESS_END                    = 100102
const GLU_END                         = 100102
const GLU_TESS_ERROR                  = 100103
const GLU_TESS_EDGE_FLAG              = 100104
const GLU_EDGE_FLAG                   = 100104
const GLU_TESS_COMBINE                = 100105
const GLU_TESS_BEGIN_DATA             = 100106
const GLU_TESS_VERTEX_DATA            = 100107
const GLU_TESS_END_DATA               = 100108
const GLU_TESS_ERROR_DATA             = 100109
const GLU_TESS_EDGE_FLAG_DATA         = 100110
const GLU_TESS_COMBINE_DATA           = 100111
const GLU_CW                          = 100120
const GLU_CCW                         = 100121
const GLU_INTERIOR                    = 100122
const GLU_EXTERIOR                    = 100123
const GLU_UNKNOWN                     = 100124
const GLU_TESS_WINDING_RULE           = 100140
const GLU_TESS_BOUNDARY_ONLY          = 100141
const GLU_TESS_TOLERANCE              = 100142
const GLU_TESS_ERROR1                 = 100151
const GLU_TESS_ERROR2                 = 100152
const GLU_TESS_ERROR3                 = 100153
const GLU_TESS_ERROR4                 = 100154
const GLU_TESS_ERROR5                 = 100155
const GLU_TESS_ERROR6                 = 100156
const GLU_TESS_ERROR7                 = 100157
const GLU_TESS_ERROR8                 = 100158
const GLU_TESS_MISSING_BEGIN_POLYGON  = 100151
const GLU_TESS_MISSING_BEGIN_CONTOUR  = 100152
const GLU_TESS_MISSING_END_POLYGON    = 100153
const GLU_TESS_MISSING_END_CONTOUR    = 100154
const GLU_TESS_COORD_TOO_LARGE        = 100155
const GLU_TESS_NEED_COMBINE_CALLBACK  = 100156
const GLU_TESS_WINDING_ODD            = 100130
const GLU_TESS_WINDING_NONZERO        = 100131
const GLU_TESS_WINDING_POSITIVE       = 100132
const GLU_TESS_WINDING_NEGATIVE       = 100133
const GLU_TESS_WINDING_ABS_GEQ_TWO    = 100134
const GLU_TESS_MAX_COORD              = 1.0e150
export GLU_EXT_object_space_tess
export GLU_EXT_nurbs_tessellator
export GLU_FALSE
export GLU_TRUE
export GLU_VERSION_1_1
export GLU_VERSION_1_2
export GLU_VERSION_1_3
export GLU_VERSION
export GLU_EXTENSIONS
export GLU_INVALID_ENUM
export GLU_INVALID_VALUE
export GLU_OUT_OF_MEMORY
export GLU_INCOMPATIBLE_GL_VERSION
export GLU_INVALID_OPERATION
export GLU_OUTLINE_POLYGON
export GLU_OUTLINE_PATCH
export GLU_NURBS_ERROR
export GLU_ERROR
export GLU_NURBS_BEGIN
export GLU_NURBS_BEGIN_EXT
export GLU_NURBS_VERTEX
export GLU_NURBS_VERTEX_EXT
export GLU_NURBS_NORMAL
export GLU_NURBS_NORMAL_EXT
export GLU_NURBS_COLOR
export GLU_NURBS_COLOR_EXT
export GLU_NURBS_TEXTURE_COORD
export GLU_NURBS_TEX_COORD_EXT
export GLU_NURBS_END
export GLU_NURBS_END_EXT
export GLU_NURBS_BEGIN_DATA
export GLU_NURBS_BEGIN_DATA_EXT
export GLU_NURBS_VERTEX_DATA
export GLU_NURBS_VERTEX_DATA_EXT
export GLU_NURBS_NORMAL_DATA
export GLU_NURBS_NORMAL_DATA_EXT
export GLU_NURBS_COLOR_DATA
export GLU_NURBS_COLOR_DATA_EXT
export GLU_NURBS_TEXTURE_COORD_DATA
export GLU_NURBS_TEX_COORD_DATA_EXT
export GLU_NURBS_END_DATA
export GLU_NURBS_END_DATA_EXT
export GLU_NURBS_ERROR1
export GLU_NURBS_ERROR2
export GLU_NURBS_ERROR3
export GLU_NURBS_ERROR4
export GLU_NURBS_ERROR5
export GLU_NURBS_ERROR6
export GLU_NURBS_ERROR7
export GLU_NURBS_ERROR8
export GLU_NURBS_ERROR9
export GLU_NURBS_ERROR10
export GLU_NURBS_ERROR11
export GLU_NURBS_ERROR12
export GLU_NURBS_ERROR13
export GLU_NURBS_ERROR14
export GLU_NURBS_ERROR15
export GLU_NURBS_ERROR16
export GLU_NURBS_ERROR17
export GLU_NURBS_ERROR18
export GLU_NURBS_ERROR19
export GLU_NURBS_ERROR20
export GLU_NURBS_ERROR21
export GLU_NURBS_ERROR22
export GLU_NURBS_ERROR23
export GLU_NURBS_ERROR24
export GLU_NURBS_ERROR25
export GLU_NURBS_ERROR26
export GLU_NURBS_ERROR27
export GLU_NURBS_ERROR28
export GLU_NURBS_ERROR29
export GLU_NURBS_ERROR30
export GLU_NURBS_ERROR31
export GLU_NURBS_ERROR32
export GLU_NURBS_ERROR33
export GLU_NURBS_ERROR34
export GLU_NURBS_ERROR35
export GLU_NURBS_ERROR36
export GLU_NURBS_ERROR37
export GLU_AUTO_LOAD_MATRIX
export GLU_CULLING
export GLU_SAMPLING_TOLERANCE
export GLU_DISPLAY_MODE
export GLU_PARAMETRIC_TOLERANCE
export GLU_SAMPLING_METHOD
export GLU_U_STEP
export GLU_V_STEP
export GLU_NURBS_MODE
export GLU_NURBS_MODE_EXT
export GLU_NURBS_TESSELLATOR
export GLU_NURBS_TESSELLATOR_EXT
export GLU_NURBS_RENDERER
export GLU_NURBS_RENDERER_EXT
export GLU_OBJECT_PARAMETRIC_ERROR
export GLU_OBJECT_PARAMETRIC_ERROR_EXT
export GLU_OBJECT_PATH_LENGTH
export GLU_OBJECT_PATH_LENGTH_EXT
export GLU_PATH_LENGTH
export GLU_PARAMETRIC_ERROR
export GLU_DOMAIN_DISTANCE
export GLU_MAP1_TRIM_2
export GLU_MAP1_TRIM_3
export GLU_POINT
export GLU_LINE
export GLU_FILL
export GLU_SILHOUETTE
export GLU_SMOOTH
export GLU_FLAT
export GLU_NONE
export GLU_OUTSIDE
export GLU_INSIDE
export GLU_TESS_BEGIN
export GLU_BEGIN
export GLU_TESS_VERTEX
export GLU_VERTEX
export GLU_TESS_END
export GLU_END
export GLU_TESS_ERROR
export GLU_TESS_EDGE_FLAG
export GLU_EDGE_FLAG
export GLU_TESS_COMBINE
export GLU_TESS_BEGIN_DATA
export GLU_TESS_VERTEX_DATA
export GLU_TESS_END_DATA
export GLU_TESS_ERROR_DATA
export GLU_TESS_EDGE_FLAG_DATA
export GLU_TESS_COMBINE_DATA
export GLU_CW
export GLU_CCW
export GLU_INTERIOR
export GLU_EXTERIOR
export GLU_UNKNOWN
export GLU_TESS_WINDING_RULE
export GLU_TESS_BOUNDARY_ONLY
export GLU_TESS_TOLERANCE
export GLU_TESS_ERROR1
export GLU_TESS_ERROR2
export GLU_TESS_ERROR3
export GLU_TESS_ERROR4
export GLU_TESS_ERROR5
export GLU_TESS_ERROR6
export GLU_TESS_ERROR7
export GLU_TESS_ERROR8
export GLU_TESS_MISSING_BEGIN_POLYGON
export GLU_TESS_MISSING_BEGIN_CONTOUR
export GLU_TESS_MISSING_END_POLYGON
export GLU_TESS_MISSING_END_CONTOUR
export GLU_TESS_COORD_TOO_LARGE
export GLU_TESS_NEED_COMBINE_CALLBACK
export GLU_TESS_WINDING_ODD
export GLU_TESS_WINDING_NONZERO
export GLU_TESS_WINDING_POSITIVE
export GLU_TESS_WINDING_NEGATIVE
export GLU_TESS_WINDING_ABS_GEQ_TWO
export GLU_TESS_MAX_COORD

end


text:
module AdversarialDriving
    using POMDPs
    using POMDPPolicies
    using POMDPSimulators
    using AutomotiveSimulator
    using AutomotiveVisualization
    using Distributions
    using Parameters
    using Random
    using Reel
    using Cairo

    # TIDM model exports
    export  PEDESTRIAN_DEF, Noise, NoisyState, NoisyPedState, BlinkerState,
            update_veh_state, noisy_scene, noisy_entity, blinker, goals, noise,
            BlinkerVehicle, NoisyPedestrian, Disturbance, BlinkerVehicleControl,
            PedestrianControl, AdversarialPedestrian, TIDM, lane_belief, laneid,
            can_have_goal, any_collides, ego_collides, end_of_road, PolicyModel,
            Blindspot, in_blindspot, RenderableBlindspot
    include("driving_models.jl")

    # T-Intersection roadway exports
    export  roadway_lane,
            road_with_crosswalk, Crosswalk, DEFAULT_CROSSWALK_LANE,
            ped_roadway, crosswalk, ped_yields_way, ped_intersection_enter_loc, ped_intersection_exit_loc, ped_goals, ped_should_blink,
            T_intersection, Tint_goal_label, Tint_signal_right,
            Tint_roadway, Tint_yields_way, Tint_intersection_enter_loc, Tint_intersection_exit_loc, Tint_goals, Tint_should_blink,
            Tint_TIDM_template, ped_TIDM_template
    include("roadways.jl")

    # MDP
    export Agent, BlinkerVehicleAgent, NoisyPedestrianAgent, id,
           AdversarialDrivingMDP, action_probability, step_scene,
           agents, adversaries, model, sut, sutid, update_adversary!, DrivingMDP
    include("mdp.jl")

    # states
    export PEDESTRIAN_ENTITY_DIM, NoisyPedestrian_to_vec, vec_to_NoisyPedestrian_fn,
           BLINKERVEHICLE_ENTITY_DIM, BLINKERVEHICLE_EXPANDED_ENTITY_DIM, BlinkerVehicle_to_vec, vec_to_BlinkerVehicle, BlinkerVehicle_to_expanded_vec
    include("states.jl")

    #actions
    export PEDESTRIAN_DISTURBANCE_DIM, PedestrianControl_to_vec, vec_to_PedestrianControl,
           BLINKERVEHICLE_DISTURBANCE_DIM, BlinkerVehicleControl_to_vec, vec_to_BlinkerVehicleControl,
           BV_ACTIONS, BV_ACTION_PROB, get_actions, get_actionindex, DiscreteActionModel,
           combine_disturbance_models, combine_continuous, combine_discrete, get_bv_actions, get_rand_bv_actions
    include("actions.jl")

    #helpers
    export decompose_scene,
           decompose_indices,
           ez_Tint_vehicle,
           get_Tint_vehicle,
           get_rand_Tint_vehicle,
           up_left,
           left_straight,
           left_turnright,
           right_straight,
           right_turnleft,
           rand_up_left,
           rand_left,
           rand_right,
           ez_ped_vehicle,
           ez_pedestrian,
           get_ped_vehicle,
           get_rand_ped_vehicle,
           get_pedestrian,
           get_rand_pedestrian,
           scenes_to_gif,
           random_IDM,
           create_actions_1BV,
           continuous_policy
    include("helpers.jl")
end

ext: jl


text:
using JuliaInterpreter
using JuliaInterpreter: Frame, @lookup
using JuliaInterpreter: finish_and_return!, evaluate_call!, step_expr!, shouldbreak,
                        do_assignment!, SSAValue, isassign, pc_expr, handle_err, get_return,
                        moduleof
using Base.Meta: isexpr
using Test, Random, SHA

function stacklength(frame)
    n = 1
    frame = frame.callee
    while frame !== nothing
        n += 1
        frame = frame.callee
    end
    return n
end

# Execute a frame using Julia's regular compiled-code dispatch for any :call expressions
runframe(frame) = Some{Any}(finish_and_return!(Compiled(), frame))

# Execute a frame using the interpreter for all :call expressions (except builtins & intrinsics)
runstack(frame) = Some{Any}(finish_and_return!(frame))

## For juliatests.jl

function read_and_parse(filename)
    src = read(filename, String)
    ex = Base.parse_input_line(src; filename=filename)
end

## For running interpreter frames under resource limitations

struct Aborted    # for signaling that some statement or test blocks were interrupted
    at::Core.LineInfoNode
end

function Aborted(frame::Frame, pc)
    src = frame.framecode.src
    lineidx = src.codelocs[pc]
    lineinfo = JuliaInterpreter.linetable(frame, lineidx; macro_caller=true)
    return Aborted(lineinfo)
end

"""
    ret, nstmtsleft = evaluate_limited!(recurse, frame, nstmts, istoplevel::Bool=true)

Run `frame` until one of:
- execution terminates normally (`ret = Some{Any}(val)`, where `val` is the returned value of `frame`)
- if `istoplevel` and a `thunk` or `method` expression is encountered (`ret = nothing`)
- more than `nstmts` have been executed (`ret = Aborted(lin)`, where `lnn` is the `LineInfoNode` of termination).
"""
function evaluate_limited!(@nospecialize(recurse), frame::Frame, nstmts::Int, istoplevel::Bool=false)
    refnstmts = Ref(nstmts)
    limexec!(s, f, istl) = limited_exec!(s, f, refnstmts, istl)
    # The following is like finish!, except we intercept :call expressions so that we can run them
    # with limexec! rather than the default finish_and_return!
    pc = frame.pc
    while nstmts > 0
        shouldbreak(frame, pc) && return BreakpointRef(frame.framecode, pc), refnstmts[]
        stmt = pc_expr(frame, pc)
        if isa(stmt, Expr)
            if stmt.head == :call && !isa(recurse, Compiled)
                refnstmts[] = nstmts
                try
                    rhs = evaluate_call!(limexec!, frame, stmt)
                    isa(rhs, Aborted) && return rhs, refnstmts[]
                    lhs = SSAValue(pc)
                    do_assignment!(frame, lhs, rhs)
                    new_pc = pc + 1
                catch err
                    new_pc = handle_err(recurse, frame, err)
                end
                nstmts = refnstmts[]
            elseif stmt.head == :(=) && isexpr(stmt.args[2], :call) && !isa(recurse, Compiled)
                refnstmts[] = nstmts
                try
                    rhs = evaluate_call!(limexec!, frame, stmt.args[2])
                    isa(rhs, Aborted) && return rhs, refnstmts[]
                    do_assignment!(frame, stmt.args[1], rhs)
                    new_pc = pc + 1
                catch err
                    new_pc = handle_err(recurse, frame, err)
                end
                nstmts = refnstmts[]
            elseif istoplevel && stmt.head == :thunk
                code = stmt.args[1]
                if length(code.code) == 1 && JuliaInterpreter.is_return(code.code[end]) && isexpr(code.code[end].args[1], :method)
                    # Julia 1.2+ puts a :thunk before the start of each method
                    new_pc = pc + 1
                else
                    refnstmts[] = nstmts
                    newframe = Frame(moduleof(frame), stmt)
                    if isa(recurse, Compiled)
                        finish!(recurse, newframe, true)
                    else
                        newframe.caller = frame
                        frame.callee = newframe
                        ret = limited_exec!(recurse, newframe, refnstmts, istoplevel)
                        isa(ret, Aborted) && return ret, refnstmts[]
                        frame.callee = nothing
                    end
                    JuliaInterpreter.recycle(newframe)
                    # Because thunks may define new methods, return to toplevel
                    frame.pc = pc + 1
                    return nothing, refnstmts[]
                end
            elseif istoplevel && stmt.head == :method && length(stmt.args) == 3
                step_expr!(recurse, frame, stmt, istoplevel)
                frame.pc = pc + 1
                return nothing, nstmts - 1
            else
                new_pc = step_expr!(recurse, frame, stmt, istoplevel)
                nstmts -= 1
            end
        else
            new_pc = step_expr!(recurse, frame, stmt, istoplevel)
            nstmts -= 1
        end
        (new_pc === nothing || isa(new_pc, BreakpointRef)) && break
        pc = frame.pc = new_pc
    end
    # Handle the return
    stmt = pc_expr(frame, pc)
    if nstmts == 0 && !JuliaInterpreter.is_return(stmt)
        ret = Aborted(frame, pc)
        return ret, nstmts
    end
    ret = get_return(frame)
    return Some{Any}(ret), nstmts
end

evaluate_limited!(@nospecialize(recurse), modex::Tuple{Module,Expr,Frame}, nstmts::Int, istoplevel::Bool=true) =
    evaluate_limited!(recurse, modex[end], nstmts, istoplevel)
evaluate_limited!(@nospecialize(recurse), modex::Tuple{Module,Expr,Expr}, nstmts::Int, istoplevel::Bool=true) =
    Some{Any}(Core.eval(modex[1], modex[3])), nstmts

evaluate_limited!(frame::Union{Frame, Tuple}, nstmts::Int, istoplevel::Bool=false) =
    evaluate_limited!(finish_and_return!, frame, nstmts, istoplevel)

function limited_exec!(@nospecialize(recurse), newframe, refnstmts, istoplevel)
    ret, nleft = evaluate_limited!(recurse, newframe, refnstmts[], istoplevel)
    refnstmts[] = nleft
    return isa(ret, Aborted) ? ret : something(ret)
end

### Functions needed on workers for running tests

function configure_test()
    # To run tests efficiently, certain methods must be run in Compiled mode,
    # in particular those that are used by the Test infrastructure
    cm = JuliaInterpreter.compiled_methods
    empty!(cm)
    JuliaInterpreter.set_compiled_methods()
    push!(cm, which(Test.eval_test, Tuple{Expr, Expr, LineNumberNode}))
    push!(cm, which(Test.get_testset, Tuple{}))
    push!(cm, which(Test.push_testset, Tuple{Test.AbstractTestSet}))
    push!(cm, which(Test.pop_testset, Tuple{}))
    for f in (Test.record, Test.finish)
        for m in methods(f)
            push!(cm, m)
        end
    end
    push!(cm, which(Random.seed!, Tuple{Union{Integer,Vector{UInt32}}}))
    push!(cm, which(copy!, Tuple{Random.MersenneTwister, Random.MersenneTwister}))
    push!(cm, which(copy, Tuple{Random.MersenneTwister}))
    push!(cm, which(Base.include, Tuple{Module, String}))
    push!(cm, which(Base.show_backtrace, Tuple{IO, Vector}))
    push!(cm, which(Base.show_backtrace, Tuple{IO, Vector{Any}}))
    # issue #101
    push!(cm, which(SHA.update!, Tuple{SHA.SHA1_CTX,Vector{UInt8}}))
end

function run_test_by_eval(test, fullpath, nstmts)
    Core.eval(Main, Expr(:toplevel, :(module JuliaTests using Test, Random end), quote
        # These must be run at top level, so we can't put this in a function
        println("Working on ", $test, "...")
        ex = read_and_parse($fullpath)
        isexpr(ex, :error) && @error "error parsing $($test): $ex"
        aborts = Aborted[]
        ts = Test.DefaultTestSet($test)
        Test.push_testset(ts)
        current_task().storage[:SOURCE_PATH] = $fullpath
        modexs = collect(ExprSplitter(JuliaTests, ex))
        for (i, modex) in enumerate(modexs)  # having the index can be useful for debugging
            nstmtsleft = $nstmts
            # mod, ex = modex
            # @show mod ex
            frame = Frame(modex)
            yield()  # allow communication between processes
            ret, nstmtsleft = evaluate_limited!(frame, nstmtsleft, true)
            if isa(ret, Aborted)
                push!(aborts, ret)
                JuliaInterpreter.finish_stack!(Compiled(), frame, true)
            end
        end
        println("Finished ", $test)
        return ts, aborts
    end))
end


text:
module DeflateTruss

using Ferrite
using StaticArrays
using NearestNeighbors
using TopOpt.TopOptProblems: StiffnessTopOptProblem, Metadata, RectilinearGrid,
    left, right, middley, middlez, find_black_and_white, find_varind
using TopOpt.TrussTopOptProblems
using TopOpt.TrussTopOptProblems: TrussProblem, get_fixities_node_set_name

export CustomPointLoadCantileverTruss

######################

function CustomPointLoadCantileverTruss(nels::NTuple{dim,Int}, sizes::NTuple{dim}, E = 1.0, ν = 0.3, force = ones(dim); k_connect=1) where {dim}
    iseven(nels[2]) && (length(nels) < 3 || iseven(nels[3])) || throw("Grid does not have an even number of elements along the y and/or z axes.")
    _T = promote_type(eltype(sizes), typeof(E), typeof(ν), eltype(force))
    if _T <: Integer
        T = Float64
    else
        T = _T
    end
    @assert length(force) == dim

    # only for the convience of getting all the node points
    rect_grid = RectilinearGrid(Val{:Linear}, nels, T.(sizes))
    node_mat = hcat(map(x -> Vector(x.x), rect_grid.grid.nodes)...)
    kdtree = KDTree(node_mat)
    if dim == 2
        # 4+1*4 -> 4+3*4 -> 4+5*4
        k_ = 4*k_connect + 4*sum(1:2:(2*k_connect-1))
    else
        k_ = 8*k_connect + 6*sum(1:9:(9*k_connect-1))
    end
    idxs, _ = knn(kdtree, node_mat, k_+1, true)
    connect_mat = zeros(Int, 2, k_*length(idxs))
    for (i, v) in enumerate(idxs)
        connect_mat[1, (i-1)*k_+1:i*k_] = ones(Int, k_)*i
        connect_mat[2, (i-1)*k_+1:i*k_] = v[2:end] # skip the point itself
    end
    truss_grid = TrussGrid(node_mat, connect_mat)

    # reference domain dimension for a line element
    ξdim = 1
    ncells = getncells(truss_grid)
    mats = [TrussFEAMaterial{T}(E, ν) for i=1:ncells]

    # * support nodeset
    for i in 1:dim
        addnodeset!(truss_grid.grid, get_fixities_node_set_name(i), x -> left(rect_grid, x));
    end

    # * load nodeset
    if dim == 2
        addnodeset!(truss_grid.grid, "force", x -> right(rect_grid, x) && middley(rect_grid, x));
    else
        addnodeset!(truss_grid.grid, "force", x -> right(rect_grid, x) && middley(rect_grid, x)
            && middlez(rect_grid, x));
    end

    # * Create displacement field u
    geom_order = 1
    dh = DofHandler(truss_grid.grid)
    ip = Lagrange{ξdim, RefCube, geom_order}()
    push!(dh, :u, dim, ip)
    close!(dh)

    ch = ConstraintHandler(dh)
    for i in 1:dim
        dbc = Dirichlet(
            :u,
            getnodeset(truss_grid.grid, get_fixities_node_set_name(i)),
            (x, t) -> zeros(T, 1),
            [i],
        )
        add!(ch, dbc)
    end
    close!(ch)
    # update the DBC to current time
    t = T(0)
    update!(ch, t)

    metadata = Metadata(dh)

    loadset = getnodeset(truss_grid.grid, "force")
    ploads = Dict{Int, SVector{dim, T}}()
    for node_id in loadset
        ploads[node_id] = SVector{dim,T}(dim == 2 ? force : force)
    end

    black, white = find_black_and_white(dh)
    varind = find_varind(black, white)

    return TrussProblem(truss_grid, mats, ch, ploads, black, white, varind, metadata)
end

end # end module


text:
"""
$(TYPEDEF)
"""
mutable struct Point2Point2Range{D <: IIF.SamplableBelief} <: IncrementalInference.FunctorPairwiseMinimize
  Z::D
  Point2Point2Range{D}() where {D} = new{D}()
  Point2Point2Range{D}(d::D) where {D <: IIF.SamplableBelief} = new{D}(d)
end
Point2Point2Range(d::D) where {D <: IIF.SamplableBelief} = Point2Point2Range{D}(d)
function getSample(pp2::Point2Point2Range{T}, N::Int=1) where {T <: IIF.SamplableBelief}
  return (reshape(rand(pp2.Z,N),1,N),  2*pi*rand(N))
end
function (pp2r::Point2Point2Range{T})(
            res::Array{Float64},
            userdata::FactorMetadata,
            idx::Int,
            meas::Tuple,
            xi::Array{Float64,2},
            lm::Array{Float64,2} ) where {T <: IIF.SamplableBelief}
  #
  z = meas[1][1,idx]
  XX = lm[1,idx] - (z*cos(meas[2][idx]) + xi[1,idx])
  YY = lm[2,idx] - (z*sin(meas[2][idx]) + xi[2,idx])
  res[1] = XX^2 + YY^2
  # @show "fnc", pointer(res), res
  res[1]
end
# import RoME: Point2Point2Range

passTypeThrough(d::FunctionNodeData{Point2Point2Range}) = d

"""
$(TYPEDEF)
"""
mutable struct PackedPoint2Point2Range  <: IncrementalInference.PackedInferenceType
  str::String
  PackedPoint2Point2Range() = new()
  PackedPoint2Point2Range(s::AS) where {AS <: AbstractString} = new(s)
end
function convert(::Type{PackedPoint2Point2Range}, d::Point2Point2Range)
  return PackedPoint2Point2Range(string(d.Z))
end
function convert(::Type{Point2Point2Range}, d::PackedPoint2Point2Range)
  return Point2Point2Range(extractdistribution(d.str))
end

"""
    $TYPEDEF

Range only measurement from Pose2 to Point2 variable.
"""
mutable struct Pose2Point2Range{T} <: IncrementalInference.FunctorPairwiseMinimize
  Z::T
  partial::Tuple{Int,Int}
  Pose2Point2Range{T}() where T = new()
  Pose2Point2Range{T}(Z::T) where {T <: IIF.SamplableBelief} = new{T}(Z, (1,2))
end
Pose2Point2Range(Z::T) where {T <: IIF.SamplableBelief} = Pose2Point2Range{T}(Z)

function getSample(pp2::Pose2Point2Range, N::Int=1)
  return (reshape(rand(pp2.Z,N),1,N) ,  2*pi*rand(N))
end
function (pp2r::Pose2Point2Range)(res::Array{Float64},
                                  userdata,
                                  idx::Int,
                                  meas::Tuple{Array{Float64,2}, Array{Float64,1}}, # from getSample
                                  xi::Array{Float64,2},
                                  lm::Array{Float64,2}  )
  #
  # DONE in IIF -- still need to add multi-hypotheses support here
  # this is the noisy range
  z = meas[1][1,idx]
  XX = lm[1,idx] - (z*cos(meas[2][idx]) + xi[1,idx])
  YY = lm[2,idx] - (z*sin(meas[2][idx]) + xi[2,idx])
  res[1] = XX^2 + YY^2
  return res[1]
end

mutable struct PackedPose2Point2Range  <: IncrementalInference.PackedInferenceType
  str::String
  PackedPose2Point2Range() = new()
  PackedPose2Point2Range(s::AS) where {AS <: AbstractString} = new(s)
end
function convert(::Type{PackedPose2Point2Range}, d::Pose2Point2Range)
  return PackedPose2Point2Range(string(d.Z))
end
function convert(::Type{Pose2Point2Range}, d::PackedPose2Point2Range)
  return Pose2Point2Range(extractdistribution(d.str))
end


text:
include("solver_preamble.jl")
import GLPK
factory = GLPK.Optimizer
config = MOI.Test.TestConfig(atol=1e-5, rtol=1e-5, query=false)
@testset "Linear" begin
    # With `dsos_horn_test`, the termination status is
    # `INFEASIBLE_OR_UNBOUNDED` instead of `INFEASIBLE`.
    Tests.linear_test(factory, config, ["dsos_horn"])
end


text:
include("double_learner.jl")
include("monte_carlo_learner.jl")
include("temporal_difference_learner.jl")
include("gradient_bandit_learner.jl")
include("patch.jl")


text:
"""
    modularity(g, c, distmx=weights(g), γ=1.0)

Return a value representing Newman's modularity `Q` for the undirected and
directed graph `g` given the partitioning vector `c`. This method also supports
weighted graphs if the distance matrix is provided.

Modularity ``Q`` for undirected graph:
```math
Q = \\frac{1}{2m} \\sum_{c} \\left( e_{c} - \\gamma \\frac{K_c^2}{2m} \\right)
```

Modularity ``Q`` for directed graph:
```math
Q = \\frac{1}{m} \\sum_{c} \\left( e_{c} - \\gamma \\frac{K_c^{in} K_c^{out}}{m} \\right)
```

where:
- ``m``: total number of edges in the network
- ``e_c``: number of edges in community ``c``
- ``K_c``: sum of the degrees of the nodes in community ``c`` or the
  sum of the weighted degree of the nodes in community ``c`` when the graph is
  weighted. ``K_c^{in}`` sum of the in-degrees of the nodes in community ``c``.

### Optional Arguments
- `distmx=weights(g)`: distance matrix for weighted graphs
- `γ=1.0`: where `γ > 0` is a resolution parameter. When the modularity is used
  to find communities structure in networks (i.e with [Louvain's method for
  community detection](https://en.wikipedia.org/wiki/Louvain_Modularity)),
  higher resolutions lead to more communities, while lower resolutions lead to
  fewer communities. Where `γ=1.0` it lead to the traditional definition of
  the modularity.

### References
- M. E. J. Newman and M. Girvan. "Finding and evaluating community structure in networks".
  Phys. Rev. E 69, 026113 (2004). [(arXiv)](https://arxiv.org/abs/cond-mat/0308217)
- J. Reichardt and S. Bornholdt. "Statistical mechanics of community detection".
  Phys. Rev. E 74, 016110 (2006). [(arXiv)](https://arxiv.org/abs/cond-mat/0603718)
- E. A. Leicht and M. E. J. Newman. "Community structure in directed networks".
  Physical Review Letter, 100:118703, (2008). [(arXiv)](https://arxiv.org/pdf/0709.4500.pdf)

# Examples
```jldoctest
julia> using LightGraphs

julia> barbell = blockdiag(complete_graph(3), complete_graph(3));

julia> add_edge!(barbell, 1, 4);

julia> modularity(barbell, [1, 1, 1, 2, 2, 2])
0.35714285714285715

julia> modularity(barbell, [1, 1, 1, 2, 2, 2], γ=0.5)
0.6071428571428571

julia> using SimpleWeightedGraphs

julia> triangle = SimpleWeightedGraph(3);

julia> add_edge!(triangle, 1, 2, 1);

julia> add_edge!(triangle, 2, 3, 1);

julia> add_edge!(triangle, 3, 1, 1);

julia> barbell = blockdiag(triangle, triangle);

julia> add_edge!(barbell, 1, 4, 5); # this edge has a weight of 5

julia> modularity(barbell, [1, 1, 1, 2, 2, 2])
0.045454545454545456
```
"""
function modularity(
    g::AbstractGraph,
    c::AbstractVector{<:Integer};
    distmx::AbstractArray{<:Number}=weights(g),
    γ=1.0
    )

    m = sum([distmx[src(e), dst(e)] for e in edges(g)])
    m = is_directed(g) ? m : 2 * m

    m == 0 && return 0.
    nc = maximum(c)
    kin = zeros(Float32, nc)
    kout = zeros(Float32, nc)
    Q = 0.0
    for u in vertices(g)
        for v in neighbors(g, u)
            c1 = c[u]
            c2 = c[v]
            if c1 == c2
                Q += distmx[u,v]
            end
            kout[c1] += distmx[u,v]
            kin[c2] += distmx[u,v]
        end
    end
    Q = Q * m
    @inbounds for i = 1:nc
        Q -= γ * kin[i] * kout[i]
    end
    return Q / m^2
end


text: println("Model Build Successful!")


text:
#EXERCISE MADE BY ERIK MARTÍN GARZÓN & DAVID BOCES OBIS
# the packages we need
using Gtk, Graphics, Logging, Printf

include("proj2.jl")

# the main window
win = GtkWindow("SO(3)")

function init_canvas(h,w)
    # create the drawing canvas
    c = @GtkCanvas(h,w)

    # create the initial drawing inside a try/catch loop
    @guarded draw(c) do widget
        #  draw the background with the canvas drawing context
        # the code for this comes later
        draw_the_canvas(c)
    end

    show(c)
    return c
end

# make the canvas
the_canvas = init_canvas(750,750)

# --------- part 2 -------------
# define all the widgets

# a widget for status messages that we define at the beginning so we can use it from the callback
msg_label = GtkLabel("No message at this time")

matrix_label1 = GtkLabel("")
matrix_label2 = GtkLabel("")
matrix_label3 = GtkLabel("")

global matrixR = [0.0 0.0 0.0; 0.0 0.0 0.0; 0.0 0.0 0.0]

# defaults
default_value = Dict("phi"=>0,"v_x"=>1,"v_y"=>0,"v_z"=>0,"alpha"=>70, "q0" => 0, "q1" => 0, "q2" => 0, "q3" => 0, "m1" => 0, "m2" => 0, "m3" => 0)

# an array to store the entry boxes
entry_list = []

# an array of labels that we use to display normalized inputs,
# and which also gets modified from the callback
normalized_labels = []

#we comprobate if quaternions change, for now, is false
quatChanged = false

function find_by_name(list, name)
    for item in list
        if get_gtk_property(item, :name, String) == name
            return item
        end
    end
    @warn name, "not found in list"
    @warn "available names are"
    for item in list
        @warn get_gtk_property(item, :name, String)
    end
end

function output_normalized(label, value)
    GAccessor.text(find_by_name(normalized_labels, label), @sprintf("%3.2f", value))
end

function normalize_v()
    v_x = read_original_box("v_x")
    v_y = read_original_box("v_y")
    v_z = read_original_box("v_z")

    norm = sqrt(v_x*v_x + v_y*v_y + v_z*v_z)

    if(norm!=0)
        output_normalized("v_x_normalized", v_x / norm)
        output_normalized("v_y_normalized", v_y / norm)
        output_normalized("v_z_normalized", v_z / norm)
    else
        output_normalized("v_x_normalized", 0)
        output_normalized("v_y_normalized", 0)
        output_normalized("v_z_normalized", 0)
    end
end

function normalize_q()
    q0 = read_original_box("q0")
    q1 = read_original_box("q1")
    q2 = read_original_box("q2")
    q3 = read_original_box("q3")

    normQ = sqrt(q0*q0 + q1*q1 + q2*q2 + q3*q3)

    if(normQ == 0)
        output_normalized("q0_normalized", 0)
        output_normalized("q1_normalized", 0)
        output_normalized("q2_normalized", 0)
        output_normalized("q3_normalized", 0)
    else
        output_normalized("q0_normalized", q0 / normQ)
        output_normalized("q1_normalized", q1 / normQ)
        output_normalized("q2_normalized", q2 / normQ)
        output_normalized("q3_normalized", q3 / normQ)

    end

end

function normalize_alpha()
    output_normalized("alpha_normalized", read_original_box("alpha"))
end

function normalize_phi()
    output_normalized("phi_normalized", read_original_box("phi"))
end

function matrix_box1()
    GAccessor.text(matrix_label1, string(matrixR[1,:]))
end

function matrix_box2()
    GAccessor.text(matrix_label2, string(matrixR[2,:]))
end

function matrix_box3()
    GAccessor.text(matrix_label3, string(matrixR[3,:]))
end

function entry_box_callback(widget)
    # who called us?
    name = get_gtk_property(widget, :name, String)
    text = get_gtk_property(widget, :text, String)

    get_gtk_property(widget, :name, String) == "matrixR"

    if (get_gtk_property(widget, :name, String) == "q0" ||
        get_gtk_property(widget, :name, String) == "q1" ||
        get_gtk_property(widget, :name, String) == "q2" ||
        get_gtk_property(widget, :name, String) == "q3")

        global quatChanged = true

    elseif (get_gtk_property(widget, :name, String) == "alpha" ||
            get_gtk_property(widget, :name, String) == "v_x" ||
            get_gtk_property(widget, :name, String) == "v_y" ||
            get_gtk_property(widget, :name, String) == "v_z")
        global quatChanged = false
    end

    # checking that we tell user alpha has limits at 0 and 180
    if (get_gtk_property(widget, :name, String) == "alpha" && read_original_box("alpha") < 0)
        GAccessor.text(msg_label, "alpha" * " changed to 0")

    elseif (get_gtk_property(widget, :name, String) == "alpha" && read_original_box("alpha") > 180)
        GAccessor.text(msg_label, "alpha" * " changed to 180" )

    else
        GAccessor.text(msg_label, name * " changed to " * text)
    end

    # change the correct normalized output
    if name[1] == 'v'
        normalize_v()
        matrix_box1()
        matrix_box2()
        matrix_box3()
    elseif name[1] == 'a'
        normalize_alpha()
        matrix_box1()
        matrix_box2()
        matrix_box3()
    elseif name[1] == 'p'
        normalize_phi()
        matrix_box1()
        matrix_box2()
        matrix_box3()
    elseif name[1] == 'q'
        normalize_q()
        matrix_box1()
        matrix_box2()
        matrix_box3()
    end

    # actually draw the changes
    draw_the_canvas(the_canvas)
    reveal(the_canvas)
end

function entry_box(label_string)
    # set up the entry
    entry = GtkEntry()
    set_gtk_property!(entry,:width_chars, 5)
    set_gtk_property!(entry,:max_length, 5)
    set_gtk_property!(entry,:name, label_string)

    default_text = string(default_value[label_string])
    GAccessor.text(entry, default_text)
    push!(entry_list, entry)

    # make it communicate changes
    signal_connect(entry_box_callback, entry, "changed")

    # set up the label and normalized output
    label = GtkLabel(label_string)
    normalized_output = GtkLabel(default_text)
    set_gtk_property!(normalized_output, :name, label_string * "_normalized")

    # make and return the containing box
    hbox = GtkButtonBox(:h)
    push!(hbox, label)
    push!(hbox, entry)
    push!(hbox, normalized_output)

    # export the normalized output for further use
    push!(normalized_labels, normalized_output)

    return hbox
end

function bold_label(label_string)
    label = GtkLabel("")
    GAccessor.markup(label, """<b>""" * label_string * """</b>""")
    return label
end

function phi_box()
    vbox = GtkBox(:v)
    push!(vbox, bold_label("Coordinate rotation"))
    push!(vbox, entry_box("phi"))
    return vbox
end

function vector_angle_box()
    vbox = GtkBox(:v)

    push!(vbox, bold_label("Axis"))

    for label in ["v_x", "v_y", "v_z"]
        push!(vbox, entry_box(label))
    end

    push!(vbox, bold_label("Angle"))

    push!(vbox, entry_box("alpha"))

    return vbox
end

function quaternion_box()
    vbox = GtkBox(:v)

    push!(vbox, bold_label("Quaternion"))
    push!(vbox, entry_box("q0"))
    push!(vbox, entry_box("q1"))
    push!(vbox, entry_box("q2"))
    push!(vbox, entry_box("q3"))
end

# Now put everything into the window,
# including the canvas
function init_window(win, canvas)

    # make a vertically stacked box for the data entry widgets
    control_box = GtkBox(:v)
    push!(control_box, phi_box())
    push!(control_box, GtkLabel(""))
    push!(control_box, vector_angle_box())
    push!(control_box, GtkLabel(""))
    push!(control_box, quaternion_box())
    push!(control_box, GtkLabel(""))
    push!(control_box, bold_label("Matrix"))
    push!(control_box, matrix_box1())
    push!(control_box, GtkLabel(""))
    push!(control_box, matrix_box2())
    push!(control_box, GtkLabel(""))
    push!(control_box, matrix_box3())
    push!(control_box, GtkLabel(""))
    push!(control_box, msg_label)

    # make another box for the drawing canvas
    canvas_box = GtkBox(:v)
    push!(canvas_box, canvas)

    # make a containing box that will stack the widgets and the canvas side by side
    global_box = GtkBox(:h)
    push!(global_box, control_box)
    push!(global_box, GtkLabel("   ")) # a very basic separator
    push!(global_box, canvas_box)

    # put it all inside the window
    push!(win, global_box)
end

# --------- part 3 -------------
# now we make the canvas interactive

function read_box(name, from_which_list, what)
    the_box = find_by_name(from_which_list, name)
    result = parse(Float64, get_gtk_property(the_box, what, String))
    return result
end

function read_original_box(name)
    return read_box(name, entry_list, :text)
end

function read_normalized_label(name)
    return read_box(name, normalized_labels, :label)
end

# the background drawing
function draw_the_canvas(canvas)
    h   = height(canvas)
    w   =  width(canvas)
    ctx =  getgc(canvas)

    # clear the canvas
    rectangle(ctx, 0, 0, w, h)
    set_source_rgb(ctx, 1, 1, 1)
    fill(ctx)

    # Paint red rectangle
    # rectangle(ctx, 3w/4, 0, w/4, h/4)
    # set_source_rgb(ctx, 1, 0, 0)
    # fill(ctx)

    # Paint blue rectangle
    # rectangle(ctx, 0, 3h/4, w/4, h/4)
    # set_source_rgb(ctx, 0, 0, 1)
    # fill(ctx)

    # Paint an example line
    # set_line_width(ctx, 5)
    # set_source_rgb(ctx, 1, 0, 1)
    # move_to(ctx, 10, 10)
    # line_to(ctx, 100, 50)
    # stroke(ctx)

    # read some normalized boxes and draw a line
    phi = 5  * read_normalized_label("phi_normalized")
    alpha = read_original_box("alpha") * 0.05
    v_x = 25 * read_normalized_label("v_x_normalized")
    v_y = 25 * read_normalized_label("v_y_normalized")
    v_z = 25 * read_normalized_label("v_z_normalized")

    q0 = read_normalized_label("q0_normalized")
    q1 = read_normalized_label("q1_normalized")
    q2 = read_normalized_label("q2_normalized")
    q3 = read_normalized_label("q3_normalized")

    q = quat(q0, q1, q2, q3)

    if(quatChanged == true)
        angle_axis = quat_to_axis_angle(q)

        alpha = angle_axis[1] * 0.05
        v_x = 25 * angle_axis[2][1]
        v_y = 25 * angle_axis[2][2]
        v_z = 25 * angle_axis[2][3]

        # println(angle_axis[1])
        # println(angle_axis[2])
    else
        qConverted = axis_angle_to_quat([v_x;v_y;v_z], alpha)

        q0 = qConverted.s
        q1 = qConverted.v1
        q2 = qConverted.v2
        q3 = qConverted.v3
    end

    global matrixR = axis_angle_to_mat([v_x;v_y;v_z], alpha)

    # println(alpha," ", q0," ", q1," ", q2," ", q3)

    if(alpha <= 0)
        alpha = 0
    end

    if(alpha >= 180*0.05)
        alpha = 180*0.05
    end

    #AXIS R3
    X = [1;0;0]
    Y = [0;1;0]
    Z = [0;0;1]

    v = [v_x; v_y; v_z]

    #WE ROTATE AXIS IN R3
    Xr = rotate_phi_z(phi, Z, X)
    Yr = rotate_phi_z(phi, Z, Y)
    Zr = rotate_phi_z(phi, Z, Z)

    vr = rotate_phi_z(phi, Z, v)

    # println("Alpha ", alpha*20)
    # println("Vr*alpha ", vr*alpha)
    # println("Xr ", Xr)

    Xvr = rotate_phi_z(alpha*20, vr*alpha, Xr)
    Yvr = rotate_phi_z(alpha*20, vr*alpha, Yr)
    Zvr = rotate_phi_z(alpha*20, vr*alpha, Zr)

    # println("Xvr ", Xvr)
    # println("Yvr ", Yvr)
    # println("Zvr ", Zvr)

    #AXIS IN R3 SCALED -VECTOR-
    Xv = scale_and_translation(Xvr*100, 0.25, vr*alpha)
    Yv = scale_and_translation(Yvr*100, 0.25, vr*alpha)
    Zv = scale_and_translation(Zvr*100, 0.25, vr*alpha)

    #WE PROJECT THE ROTATED AXIS IN R3 TO R2
    Xr2 = to_2d(Xr)
    Yr2 = to_2d(Yr)
    Zr2 = to_2d(Zr)

    #WE PROJECT THE ROTATED AXIS IN R3 TO R2
    Xvr2 = to_2d(Xv)
    Yvr2 = to_2d(Yv)
    Zvr2 = to_2d(Zv)

    #WE PROJECT THE ROTATED AXIS IN R3 TO R2
    X2 = to_2d(X)
    Y2 = to_2d(Y)
    Z2 = to_2d(Z)

    vr2 = to_2d(vr)

    #println("Xv ", Xv, "\nYv ", Yv, "\nZv ", Zv)

    #DRAW X AXIS (RED)
    set_line_width(ctx, 2)
    set_source_rgb(ctx, 1, 0, 0)
    move_to(ctx, 350, 350)
    line_to(ctx, 350 + Xr2[1]*100, 350 - Xr2[2]*100)
    stroke(ctx)

    #circle at the end of the axis
    circle(ctx, 350 + Xr2[1]*100, 350 - Xr2[2]*100, 5)
    set_source_rgb(ctx, 1, 0, 0)
    fill(ctx)

    #DRAW X AXIS SCALED -VECTOR- (RED)
    set_line_width(ctx, 2)
    set_source_rgb(ctx, 1, 0, 0)
    move_to(ctx, 350 + vr2[1]*alpha, 350 - vr2[2]*alpha)
    line_to(ctx, 350 + vr2[1] + Xvr2[1], 350 - vr2[2] - Xvr2[2])
    stroke(ctx)

    #circle at the end of the axis
    circle(ctx, 350 + vr2[1] + Xvr2[1], 350 - vr2[2] - Xvr2[2], 2)
    set_source_rgb(ctx, 1, 0, 0)
    fill(ctx)

    #DRAW Y AXIS (GREEN)
    set_line_width(ctx, 2)
    set_source_rgb(ctx, 0, 1, 0)
    move_to(ctx, 350, 350)
    line_to(ctx, 350 + Yr2[1]*100, 350 - Yr2[2]*100)
    stroke(ctx)

    #circle at the end of the axis
    circle(ctx, 350 + Yr2[1]*100, 350 - Yr2[2]*100, 5)
    set_source_rgb(ctx, 0, 1, 0)
    fill(ctx)

    #DRAW Y AXIS SCALED -VECTOR- (GREEN)
    set_line_width(ctx, 2)
    set_source_rgb(ctx, 0, 1, 0)
    move_to(ctx, 350 + vr2[1]*alpha, 350 - vr2[2]*alpha)
    line_to(ctx, 350 + vr2[1] + Yvr2[1], 350 - vr2[2] - Yvr2[2])
    stroke(ctx)

    #circle at the end of the axis
    circle(ctx, 350 + vr2[1] + Yvr2[1], 350 - vr2[2] - Yvr2[2], 2)
    set_source_rgb(ctx, 0, 1, 0)
    fill(ctx)

    #DRAW Z AXIS (BLUE)
    set_line_width(ctx, 2)
    set_source_rgb(ctx, 0, 0, 1)
    move_to(ctx, 350, 350)
    line_to(ctx, 350 + Zr2[1]*100, 350 - Zr2[2]*100)
    stroke(ctx)

    #circles at the end of the axis
    circle(ctx, 350 + Zr2[1]*100, 350 - Zr2[2]*100, 5)
    set_source_rgb(ctx, 0, 0, 1)
    fill(ctx)

    #DRAW Z AXIS SCALED -VECTOR- (BLUE)
    set_line_width(ctx, 2)
    set_source_rgb(ctx, 0, 0, 1)
    move_to(ctx, 350 + vr2[1]*alpha, 350 - vr2[2]*alpha)
    line_to(ctx, 350 + vr2[1] + Zvr2[1], 350 - vr2[2] - Zvr2[2])
    stroke(ctx)

    #circle at the end of the axis
    circle(ctx, 350 + vr2[1] + Zvr2[1], 350 - vr2[2] - Zvr2[2], 2)
    set_source_rgb(ctx, 0, 0, 1)
    fill(ctx)

    #DRAW VECTOR (BLACK)
    set_line_width(ctx, 2)
    set_source_rgb(ctx, 0, 0, 0)

    if (350+vr2[1]*alpha>1000 ||350 - vr2[2]*alpha>1000 )
        move_to(ctx, 350, 350)
        line_to(ctx, 350, 350)
    else
        move_to(ctx, 350, 350)
        line_to(ctx, 350 + vr2[1]*alpha, 350 - vr2[2]*alpha)
    end
    stroke(ctx)

    #circles at the end of the axis
    circle(ctx, 350 + vr2[1]*alpha, 350 - vr2[2]*alpha, 3)
    set_source_rgb(ctx, 0, 0, 0)
    fill(ctx)

end

# -------- initialize everything ---------

# prepare and show the initial widgets
init_window(win, the_canvas)
showall(win)


text:
function finalize_estimation(obj::KalmanODModel, arg0::ObservedMeasurement, arg1::ProcessEstimate)
    return jcall(obj, "finalizeEstimation", void, (ObservedMeasurement, ProcessEstimate), arg0, arg1)
end

function get_estimate(obj::KalmanODModel)
    return jcall(obj, "getEstimate", ProcessEstimate, ())
end

function get_estimated_propagators(obj::KalmanODModel)
    return jcall(obj, "getEstimatedPropagators", Vector{AbstractIntegratedPropagator}, ())
end

function get_evolution(obj::NonLinearProcess, arg0::jdouble, arg1::RealVector, arg2::Measurement)
    return jcall(obj, "getEvolution", NonLinearEvolution, (jdouble, RealVector, Measurement), arg0, arg1, arg2)
end

function get_innovation(obj::NonLinearProcess, arg0::Measurement, arg1::NonLinearEvolution, arg2::RealMatrix)
    return jcall(obj, "getInnovation", RealVector, (Measurement, NonLinearEvolution, RealMatrix), arg0, arg1, arg2)
end


text:
import LPDM: bounds, best_lb_action, best_ub_action
using LightDarkLPDM

mutable struct LDBounds2d{S,A,O}
    lb_    ::Float64
    ub_    ::Float64
    best_lb_action_::Vec2
    best_ub_action_::Vec2
    rng::RNGVector # dummy RNG, not used

    function LDBounds2d{S,A,O}(::POMDP{S,A,O}) where {S,A,O}
        this = new{S,A,O}()
        this.lb_ = +Inf
        this.ub_ = -Inf
        this.best_lb_action_ = Vec2(NaN,NaN)
        this.best_ub_action_ = Vec2(NaN,NaN)
        this.rng = RNGVector(1, UInt32(42))
        return this
    end
end

function LPDM.bounds(b::LDBounds2d{S,A,O},
                     pomdp::AbstractLD2,
                     particles::Vector{LPDMParticle{S}},
                     config::LPDMConfig) where {S,A,O}

    # reset every time bounds are about to be recomputed
    b.lb_ = +Inf
    b.ub_ = -Inf
    b.best_lb_action_ = Vec2(NaN,NaN)
    b.best_ub_action_ = Vec2(NaN,NaN)

    tmp_lb = 0.0
    tmp_ub = 0.0
    tmp_lb_action = Vec2(0.0,0.0)
    tmp_ub_action = Vec2(0.0,0.0)

    for p in particles
        tmp_ub, tmp_ub_action = upper_bound(pomdp,p)
        tmp_lb, tmp_lb_action = lower_bound(pomdp,p) #NOTE: Option 1
        # println("s=$(p.state), lb=$tmp_lb, ub=$tmp_ub, lba=$tmp_lb_action, uba=$tmp_ub_action")

        #NOTE: Option 2
        # if p.state[1] > pomdp.min_noise_loc # assume min_noise_loc > 0
        #     # both will have to do roughly the same thing (+/- discretization differences),
        #     # so make them the same
        #     #TODO: review for lb > ub
        #     tmp_lb, tmp_lb_action = tmp_ub, tmp_ub_action
        # else
        #     tmp_lb, tmp_lb_action = lower_bound(pomdp,p)
        # end

        if tmp_lb < b.lb_
            b.lb_             = tmp_lb
            b.best_lb_action_ = tmp_lb_action
        end
        if tmp_ub > b.ub_
            b.ub_             = tmp_ub
            b.best_ub_action_ = tmp_ub_action
        end
    end

    # Sanity check
    if b.ub_ < b.lb_
        # show(particles); println("")
        # error("BOUNDS: ub=$(b.ub_) < lb=$(b.lb_), lba = $(b.best_lb_action_), uba = $(b.best_ub_action_)")
        #NOTE: let's try this to avoid discretization issues
        b.lb_ = b.ub_
    end
    return b.lb_, b.ub_
end

LPDM.best_lb_action(b::LDBounds2d) = isnan(b.best_lb_action_[1])||isnan(b.best_lb_action_[2]) ? error("best_lb_action undefined ($(b.best_lb_action_)). Call bounds() first.") : b.best_lb_action_
# LPDM.best_ub_action(b::LDBounds2d) = isnan(b.best_ub_action_) ? error("best_ub_action undefined. Call bounds() first.") : b.best_ub_action_
LPDM.best_ub_action(b::LDBounds2d) = isnan(b.best_ub_action_[1])||isnan(b.best_ub_action_[2]) ? error("best_ub_action undefined ($(b.best_ub_action_)). Call bounds() first.") : b.best_ub_action_

#TODO: these functions make sense for step-wise rewards, but not as much for quadratic rewards. May need to redo.

# computes the cost of traveling to the low noise region and only then towards target. i.e. a slow approach
   #       s__________
   #                  |
   #                  |
   #           (0,0) _|

# function lower_bound(p::AbstractLD2, particle::LPDMParticle{Vec2})
#
#     s = particle.state
#
#     r1, a1_x = coarse_move(p, s, Vec2(p.min_noise_loc,s[2]), 1)
#     r2, a1_y = coarse_move(p, Vec2(p.min_noise_loc,s[2]), Vec2(p.min_noise_loc,0.0), 2)
#     r3, a2_x = coarse_move(p, Vec2(p.min_noise_loc,0.0), Vec2(0.0,0.0), 1)
#
#     return r1 + r2 + r3, Vec2(a1_x,a1_y)
# end

# function lower_bound(p::AbstractLD2, particle::LPDMParticle{Vec2})
#
#     s = particle.state
#
#     r_right, a_right = coarse_move(p, s, Vec2(p.min_noise_loc,s[2]), 1)
#     r_down,  a_down  = coarse_move(p, Vec2(p.min_noise_loc,s[2]), Vec2(p.min_noise_loc,0.0), 2)
#     r_left,  a_left  = coarse_move(p, Vec2(p.min_noise_loc,0.0), Vec2(0.0,0.0), 1)
#
#     r_diag = r_left < r_down ? r_left : r_down
#
#     return r_right + r_diag, Vec2(a_right,0.0) #move strictly horizontally
# end

function lower_bound(p::AbstractLD2, particle::LPDMParticle{Vec2})
# println("LB===========================")
    s = particle.state

    # println("lb move 1")
    r_x1, a_x1 = move(p, s, Vec2(p.min_noise_loc,s[2]), 1)
    # println("lb move 2")
    r_y,  a_y  = move(p, Vec2(p.min_noise_loc,s[2]), Vec2(p.min_noise_loc,0.0), 2)
    # println("lb move 3")
    r_x2,  a_x2  = move(p, Vec2(p.min_noise_loc,0.0), Vec2(0.0,0.0), 1)

    # r_diag = r_x2 < r_y ? r_x2 : r_y

    if isnan(a_x1)
        # println("a_x1 = NaN")
        if isnan(a_y)
            # println("a_y = NaN")
            if isnan(a_x2)
                # println("a_x2 = NaN")
                a = Vec2(0.0,0.0) # nothing to do
            else
                a = Vec2(a_x2,0.0)
            end
        else
            a = Vec2(0.0,a_y)
        end
    else
        a = Vec2(a_x1,0.0)
    end

    return r_x1+r_y+r_x2, a
end

# function upper_bound(p::AbstractLD2, particle::LPDMParticle{Vec2})
#     s = particle.state
#
#     rx, a1_x = coarse_move(p, s, Vec2(0.0,s[2]), 1)
#     ry, a1_y = coarse_move(p, Vec2(0.0,s[2]), Vec2(0.0,0.0), 2)
#
#     # println("s=$s, rx=$rx, ry=$ry, a1_x=$a1_x, a1_y=$a1_y")
#     # error("done, for now")
#     r = rx < ry ? rx : ry                                  ## pick the smaller (worst) of the two
#     return r, Vec2(a1_x,a1_y)
# end

# computes the reward for the straight-line path to target
function upper_bound(p::AbstractLD2, particle::LPDMParticle{Vec2})
# println("UB===========================")
    s = particle.state

    # println("ub move 1")
    rx, a1_x = move(p, s, Vec2(0.0,s[2]), 1)
    # println("ub move 2")
    ry, a1_y = move(p, Vec2(0.0,s[2]), Vec2(0.0,0.0), 2)

    # println("s=$s, rx=$rx, ry=$ry, a1_x=$a1_x, a1_y=$a1_y")
    # error("done, for now")
    r = rx < ry ? rx : ry                                  ## pick the smaller (worst) of the two
    return r, Vec2(a1_x,a1_y)
end

#w1 and w2 are start and end waypoints, coord is 1 or 2 for x and y, respectively
# NOTE: this function computes a rough estimate to avoid discretization issues
function coarse_move(p::AbstractLD2, w1::Vec2, w2::Vec2, c::Int64)
    direction = w2[c] > w1[c] ? 1 : -1
    all_actions = POMDPs.actions(p, true)
    pos_actions = all_actions[all_actions .> 0]

    r = 0.0
    w = [w1[1],w1[2]]
    a = NaN
    a_dir = NaN

    a = maximum(pos_actions[pos_actions .< abs(w2[c]-w1[c])+p.term_radius]) # maximum action not exceeding Δx + p.term_radius
    # n_moves = round(Int64,(abs(w2[c]-w[c])+p.term_radius)/a
    n_moves = ceil(Int64,(abs(w2[c]-w1[c])+p.term_radius)/a)
    a_dir = direction * a

    for i in 1:n_moves
        r += reward(p, Vec2(w[1],w[2]), Vec2(a_dir,0.0)) # use current state for computing the reward; (0.0,a_dir) and (a_dir,0.0) are equivalent for reward purposes
        w[c] += a_dir # take the step
    end

    if (abs(w2[1]) <= p.term_radius) && (abs(w2[2]) <= p.term_radius)
        r += reward(p,w2,Vec2(0.0,0.0)) #termination reward
    end
    # println("move: w1=$w1, w2=$w2, n_moves = $n_moves, a_dir = $a_dir, r=$r")
    return r, a_dir #action sign depends on the direction
end

 #w1 and w2 are start and end waypoints, coord is 1 or 2 for x and y, respectively
function move(p::AbstractLD2, w1::Vec2, w2::Vec2, c::Int64)
    direction = w2[c] > w1[c] ? 1 : -1
    all_actions = POMDPs.actions(p, true)
    # I = findall(!iszero, all_actions)
    pos_actions = all_actions[all_actions .> 0]
    # nz_actions = all_actions[I]
    min_a = minimum(pos_actions)

    r = 0.0
    w = [w1[1],w1[2]]
    first_a = NaN
    a_dir = NaN

    # if (abs(w1[1]) <= p.term_radius) && (abs(w1[2]) <= p.term_radius)
    #     return reward(p,w1,Vec2(0.0,0.0)), 0.0
    # end

    if (abs(w1[c]) <= p.term_radius)
        # println("$(abs(w1[c])) <= p.term_radius")
        return reward(p,w1,Vec2(0.0,0.0)), 0.0
    end

    if abs(w2[c]-w1[c])+p.term_radius < min_a #too close to take any action
        # println("abs($(w2[c])-$(w1[c]))+$(p.term_radius) < min_a: $min_a")
        return reward(p,w1,Vec2(0.0,0.0)), 0.0
    end

    while (abs(w2[c]-w[c])+p.term_radius >= min_a) && (abs(w2[c]-w[c]) > p.term_radius)
        a = maximum(pos_actions[pos_actions .< abs(w2[c]-w[c])+p.term_radius]) # maximum action not exceeding Δx + p.term_radius
        a_dir = direction * a
        if isnan(first_a)
            first_a = a_dir # assign first action (directional)
        end
        r += reward(p, Vec2(w[1],w[2]), Vec2(a_dir,0.0)) # use current state for computing the reward; (0.0,a_dir) and (a_dir,0.0) are equivalent for reward purposes
        # println("BOUNDS: w=$w, w1=$w1, w2=$w2, min_a = $min_a, all_actions=$all_actions, pos_actions=$pos_actions, av_actions=$(all_actions[all_actions .< abs(w2[c]-w[c]) + p.term_radius]), a_dir=$a_dir,  r=$r ")
        # println("    w=($(w[1]),$(w[2])), a_dir=$a_dir, r = $r")
        w[c] += a_dir # take the step
    end

    if (abs(w2[1]) <= p.term_radius) && (abs(w2[2]) <= p.term_radius)
        # x1 < 0.7 && println("TERMINATION REWARD #2 FOR x2=$x")
        r += reward(p,w2,Vec2(0.0,0.0)) #termination reward
        # first_a = 0.0
    end

    # x1 < 0.7 && println("exiting move $x1 -> $x2")
    # if isnan(first_a) #DEBUG: remove
    #     error("move: first_a = $first_a, w1=$w1, w2=$w2")
    # end
    # println("move: w1=$w1, w2=$w2, first_a = $first_a, r=$r")
    # error("done, for now")
    return r, first_a #action sign depends on the direction
end

max_stars_repo_name: ebalaban/LightDarkLPDM

max_issues_repo_name: tomerarnon/LightDarkLPDM

max_forks_repo_name: tomerarnon/LightDarkLPDM


text:
__precompile__()
module ConvOptDL

include("data_loader.jl")
include("qp.jl")
include("svm.jl")
include("layers/dropblock.jl")
include("layers/basic_block.jl")
include("utils.jl")

end


text:
"""
    FAD(seqs; alpha = 0.01, neigh_thresh = 1.0,method = 2,err_rate=0.02, phreds = nothing)

Exploits abundance and neighborhood information to denoise sequences without clustering
on consensus calls. Primarily designed for short amplicons, high quality sequences,
and good read-per-template coverage
"""

function FAD(seqs; alpha = 0.01, neigh_thresh = 1.0,method = 2,err_rate=0.02, phreds = nothing)
    if !(method in [1,2,3])
        warn("Please pick a method in 1,2,3")
        return [],[]
    end
    #Method list:
    #1: No statistical considerations. All smaller relatives get eaten.
    #2: Use Poisson test to decide to add not-too-small relatives back in.
    #3: Use homopolymer test to skew Poisson rate for adding smaller relatives back in.

    #Need to mod this to work for .fasta files. Perhaps use the number of pairs to approximate the number of error free sequences.
    if phreds != nothing
        #seqs, phreds , _ = read_fastq(fn, seqtype=String);
        noises = [mean(phred_to_p(i)) for i in phreds];
        #Calculate the expected proportion of sequences that are error free.
        expected_zero_errors = mean([pdf(Poisson(noises[i]*length(seqs[i])),0) for i in 1:length(seqs)]);
    else
        warn("Switching to method 1, phreds missing")
        method = 1
    end

    counts = countmap(seqs);
    seq_keys = keys(counts);
    seq_freqs_all = reverse(sort([(counts[k],k,kmer_count(k,6)) for k in seq_keys]));
    seq_freqs = seq_freqs_all[[k[1]>=2 for k in seq_freqs_all]];

    if length(seq_freqs) > 0
        current_set = [seq_freqs[1]];
        for i in 2:length(seq_freqs)
            dists = [corrected_kmer_dist_full(k[3],seq_freqs[i][3]) for k in current_set]
            neighbour_set = dists.<=neigh_thresh

            #if no neighbours within radius, add to current_set
            if sum(neighbour_set)==0
                push!(current_set,seq_freqs[i])
            else
                if method == 2
                    inds = collect(1:length(current_set))[neighbour_set]
                    best = indmax(current_set[inds])
                    p_val = (1-cdf(Poisson((current_set[inds][best][1]/expected_zero_errors)*err_rate),seq_freqs[i][1]))*length(seq_freqs[i][2])
                    if p_val < alpha
                        push!(current_set,seq_freqs[i])
                    end
                end
                #Using same as above. Will add in homopolymer test.
                if method == 3
                    #Need to look through all inds. For each, calculate p-value. Take the largest p-value. If that is <thresh, then keep cluster.
                    #This will use a regular error rate if HP==false, but a much larger error rate if HP==true.

                    inds = collect(1:length(current_set))[neighbour_set]
                    best = indmax(current_set[inds])
                    p_val = (1-cdf(Poisson((current_set[inds][best][1]/expected_zero_errors)*err_rate),seq_freqs[i][1]))*length(seq_freqs[i][2])
                    if p_val < alpha
                        push!(current_set,seq_freqs[i])
                    end
                end
            end
        end

        println("Using ",100*sum([k[1] for k in current_set])/length(seqs),"% of reads. If this is lower than 5%, use RAD instead.")

        frequencies = zeros(length(current_set))
        for s in seq_freqs_all
            dists = [corrected_kmer_dist_full(k[3],s[3],k=5) for k in current_set]
            frequencies[indmin(dists)] += s[1]
        end
        return [k[2] for k in current_set],frequencies
    else #catching the condition where there are no duplicates
        warn("Failing completely. Use RAD instead.")
        return [],[]
    end
end


text:
# This currently takes the weight multiplication operator as input.
# I will probably change this to take the weight function instead.
function cholesky_jacobimatrix(w, P::OrthogonalPolynomial)
    W = Symmetric(P \ (w.(axes(P,1)) .* P)) # Compute weight multiplication via Clenshaw
    bands = CholeskyJacobiBands(W, P) # the cached array only needs to store two bands bc of symmetry
    return SymTridiagonal(bands[1,:],bands[2,:])
end

# The generated Jacobi operators are symmetric tridiagonal, so we store their data as two bands
mutable struct CholeskyJacobiBands{T} <: AbstractCachedMatrix{T}
    data::Matrix{T}
    U::UpperTriangular
    X::SymTridiagonal{T}
    datasize::Int
    array
end

# Computes the initial data for the Jacobi operator bands
function CholeskyJacobiBands(W::Symmetric{T}, P::OrthogonalPolynomial) where T
    @assert P isa Normalized
    U = cholesky(W).U
    X = jacobimatrix(P)
    dat = zeros(T,2,10)
    for k in 1:10
        dat[1,k] = (U * (X * (U \ [zeros(k-1); 1; zeros(∞)])))[k]
        dat[2,k] = (U * (X * (U \ [zeros(k); 1; zeros(∞)])))[k]
    end
    return CholeskyJacobiBands{T}(dat, U, X, 10, dat)
end

size(::CholeskyJacobiBands) = (2,ℵ₀) # Stored as two infinite bands

# Resize and filling functions for cached implementation
function resizedata!(K::CholeskyJacobiBands, nm::Integer)
    νμ = K.datasize
    if nm > νμ
        olddata = copy(K.data)
        K.data = similar(K.data, 2, nm)
        K.data[axes(olddata)...] = olddata
        inds = νμ:nm
        cache_filldata!(K, inds)
        K.datasize = nm
    end
    K
end
function cache_filldata!(J::CholeskyJacobiBands, inds::UnitRange{Int})
    for k in inds
        J.data[1,k] = (J.U * (J.X * (J.U \ [zeros(k-1); 1; zeros(∞)])))[k]
        J.data[2,k] = (J.U * (J.X * (J.U \ [zeros(k); 1; zeros(∞)])))[k]
    end
end

function getindex(K::CholeskyJacobiBands, k::Integer, j::Integer)
    resizedata!(K, max(k,j))
    K.data[k, j]
end
function getindex(K::CholeskyJacobiBands, kr::Integer, jr::UnitRange{Int})
    resizedata!(K, maximum(jr))
    K.data[kr, jr]
end
function getindex(K::CholeskyJacobiBands, I::Vararg{Int,2})
    resizedata!(K,maximum(I))
    getindex(K.data,I[1],I[2])
end


text:
module PointMassConstraintTest

using Test
using ForneyLab
using ForneyLab: outboundType, isApplicable
using ForneyLab: SPPointMassConstraintOut

#-------------
# Update rules
#-------------

@testset "SPPointMassConstraintOut" begin
    @test SPPointMassConstraintOut <: SumProductRule{PointMassConstraint}
    @test outboundType(SPPointMassConstraintOut) == Message{PointMass}
    @test isApplicable(SPPointMassConstraintOut, [Nothing])

    @test ruleSPPointMassConstraintOut(Message(Univariate, GaussianMeanVariance, m=0.0, v=1.0)) == Message(Univariate, PointMass, m=0.0)
end

#------------
# Integration
#------------

@testset "PointMassConstraint integration" begin
    fg = FactorGraph()
    @RV x ~ GaussianMeanVariance(0.0, 1.0)
    PointMassConstraint(x)

    pfz = PosteriorFactorization(fg)
    algo = messagePassingAlgorithm(x)
    algo_code = algorithmSourceCode(algo)

    @test occursin("ruleSPPointMassConstraintOut(messages[1])", algo_code)
end

end # module


text:
# This file is a part of Julia. License is MIT: https://julialang.org/license

# ImmutableArrays (arrays that implement getindex but not setindex!)

# This test file defines an array wrapper that is immutable. It can be used to
# test the action of methods on immutable arrays.

module ImmutableArrays

export ImmutableArray

"An immutable wrapper type for arrays."
struct ImmutableArray{T,N,A<:AbstractArray} <: AbstractArray{T,N}
    data::A
end

ImmutableArray(data::AbstractArray{T,N}) where {T,N} = ImmutableArray{T,N,typeof(data)}(data)

# Minimal AbstractArray interface
Base.size(A::ImmutableArray) = size(A.data)
Base.size(A::ImmutableArray, d) = size(A.data, d)
Base.getindex(A::ImmutableArray, i...) = getindex(A.data, i...)

# The immutable array remains immutable after conversion to AbstractArray
AbstractArray{T}(A::ImmutableArray) where {T} = ImmutableArray(AbstractArray{T}(A.data))
AbstractArray{T,N}(A::ImmutableArray{S,N}) where {S,T,N} = ImmutableArray(AbstractArray{T,N}(A.data))

end

max_stars_repo_name: jonas-schulze/julia

max_issues_repo_name: jonas-schulze/julia

max_forks_repo_name: jonas-schulze/julia


text:
# ---
# title: 76. Minimum Window Substring
# id: problem76
# author: Tian Jun
# date: 2020-10-31
# difficulty: Hard
# categories: Hash Table, Two Pointers, String, Sliding Window
# link: <https://leetcode.com/problems/minimum-window-substring/description/>
# hidden: true
# ---
#
# Given two strings `s` and `t`, return _the minimum window in`s` which will
# contain all the characters in `t`_. If there is no such window in `s` that
# covers all characters in `t`, return _the empty string`""`_.
#
# **Note** that If there is such a window, it is guaranteed that there will
# always be only one unique minimum window in `s`.
#
#
#
# **Example 1:**
#
#
#
#     Input: s = "ADOBECODEBANC", t = "ABC"
#     Output: "BANC"
#
#
# **Example 2:**
#
#
#
#     Input: s = "a", t = "a"
#     Output: "a"
#
#
#
#
# **Constraints:**
#
#   * `1 <= s.length, t.length <= 105`
#   * `s` and `t` consist of English letters.
#
#
#
# **Follow up:** Could you find an algorithm that runs in `O(n)` time?
#
#
## @lc code=start
using LeetCode

## add your code here:
## @lc code=end


text:
# ---
# title: 440. K-th Smallest in Lexicographical Order
# id: problem440
# author: Tian Jun
# date: 2020-10-31
# difficulty: Hard
# categories:
# link: <https://leetcode.com/problems/k-th-smallest-in-lexicographical-order/description/>
# hidden: true
# ---
#
# Given integers `n` and `k`, find the lexicographically k-th smallest integer
# in the range from `1` to `n`.
#
# Note: 1 ≤ k ≤ n ≤ 109.
#
# **Example:**
#
#
#
#     Input:
#     n: 13   k: 2
#
#     Output:
#     10
#
#     Explanation:
#     The lexicographical order is [1, 10, 11, 12, 13, 2, 3, 4, 5, 6, 7, 8, 9], so the second smallest number is 10.
#
#
#
## @lc code=start
using LeetCode

## add your code here:
## @lc code=end


text:
module Logistic

using LogicCircuits
using ..Utils

export
    LogisticΔNode,
    LogisticLeafNode,
    LogisticInnerNode,
    LogisticLiteral,
    Logistic⋀,
    Logistic⋁,
    LogisticΔ,
    LogisticΔ,
    LogisticCache,
    num_parameters_perclass,
    logistic_origin,
    class_conditional_likelihood_per_instance,
    classes

include("LogisticCircuits.jl")

end


text:
using Clang

if !isfile(joinpath(@__DIR__, "../Project.toml"))
    error("Run from scripts directory.")
end

# Download libxc sources
LIBXC_VERSION = "5.1.0"
if !isfile(joinpath(@__DIR__, "libxc/src/xc.h"))
    download("https://gitlab.com/libxc/libxc/-/archive/$LIBXC_VERSION/" *
             "libxc-$LIBXC_VERSION.tar.gz", "libxc.tar.gz")
    run(`tar xf libxc.tar.gz`)
    libxc_dir = only(d for d in readdir(".") if startswith(d, "libxc-"))
    run(`ln -sf $libxc_dir libxc`)
    rm("libxc.tar.gz")
    @assert isfile(joinpath(@__DIR__, "libxc/src/xc.h"))
end

# Generate API wrapper
LIBXC_INCLUDE = normpath(joinpath(@__DIR__, "libxc/src"))
LIBXC_HEADERS = [joinpath(LIBXC_INCLUDE, header) for header in ["xc.h"]]
OUTPUTDIR = joinpath(@__DIR__, "../src/gen")

wc = init(;headers=LIBXC_HEADERS,
          output_file=joinpath(OUTPUTDIR, "libxc.jl"),
          common_file=joinpath(OUTPUTDIR, "libxc_common.jl"),
          header_library = x -> "libxc",
          clang_includes = vcat(LIBXC_INCLUDE, CLANG_INCLUDE),
          header_wrapped = (root, current) -> root == current,
          clang_diagnostics = false,
          )

run(wc)

# Cleanup unused files
for file in ["ctypes.jl", "LibTemplate.jl"]
    rm(joinpath(OUTPUTDIR, file), force=true)
end


text:
using JuMP, EMP, ReSHOP
using Compat, Test

# Data is from  ``Data analysis recipes: Fitting a model to data'' by Hogg, Bovy & Lang
# available on https://arxiv.org/abs/1008.4686.

y = [592 401 583 402 495 173 479 504 510 416 393 442 317 311 400 337 423 334 533 344]';
x = [201 244 47 287 203 58 210 202 198 158 165 201 157 131 166 160 186 125 218 146]';
sigma_y = [61 25 38 15 21 15 27 14 30 16 14 25 52 16 34 31 42 26 16 22]';

N = length(y)
M = 1

huber_params = Dict("kappa" => 1)
l1_params = Dict()
l2_params = Dict()
elastic_net_params = Dict("lambda" => 1)
hinge_params = Dict("epsilon" => 1)
vapnik_params = Dict("epsilon" => 1)
soft_hinge_params = Dict("epsilon" => 1, "kappa" => 1)
hubnik_params = Dict("epsilon" => 1, "kappa" => 1)

all_params = Dict(
"huber" => huber_params,
"l1" => l1_params,
"l2" => l2_params,
"elastic_net" => elastic_net_params,
"hinge" => hinge_params,
"vapnik" => vapnik_params,
"soft_hinge" => soft_hinge_params,
"hubnik" => hubnik_params
)

penalty_names = keys(all_params)
ovf_formulations = ["equilibrium", "dual"]

@testset "ovf_loss_equil" begin

include(joinpath(dirname(dirname(@__FILE__)), "examples", "ovf_loss_equil.jl"))

end


text:
cd(@__DIR__) # changes the directory to the current directory, the default I guess is the HOME
using Pkg; Pkg.activate("."); Pkg.instantiate()

#= Pkg.add("NeuralPDE")
Pkg.add("Flux")
Pkg.add("ModelingToolkit")
Pkg.add("GalacticOptim")
Pkg.add("Optim")
Pkg.add("DiffEqFlux")
Pkg.add("Plots") =#

using NeuralPDE, Flux, ModelingToolkit, GalacticOptim, Optim, DiffEqFlux
import ModelingToolkit: Interval, infimum, supremum

@parameters t, x
@variables u(..)
Dxx = Differential(x)^2
Dtt = Differential(t)^2
Dt = Differential(t)

# PDE
C=1
eq = Dtt(u(t,x)) ~ C^2*Dxx(u(t,x))

# Initial & boundary conditions
bcs = [
    u(t,0) ~ 0., # for all t > 0
    u(t,1) ~ 0., # for all t > 0
    u(0,x) ~ x*(1. - x), # for all 0 < x < 1
    Dt(u(0,x)) ~ 0. # for all 0 < x < 1
]

# Space-time domains
domains = [
    t ∈ Interval(0.0, 1.0),
    x ∈ Interval(0.0, 1.0)
]

# Discretization
dx = 0.1

println("antes de chain")

# Neural Network
chain = FastChain(FastDense(2,16,Flux.σ),FastDense(16,16,Flux.σ),FastDense(16,1))
initθ = Float64.(DiffEqFlux.initial_params(chain))
discretization = PhysicsInformedNN(chain, GridTraining(dx); init_params = initθ)

@named pde_system = PDESystem(eq,bcs,domains,[t,x],[u(t,x)])
prob = discretize(pde_system,discretization)

cb = function (p,l)
    println("Current loss is: $l")
    return false
end

println("antes de optimizer")

# optimizer
opt = Optim.BFGS()
res = GalacticOptim.solve(prob,opt; cb = cb, maxiters=1200)
phi = discretization.phi

using Plots

ts,xs = [infimum(d.domain):dx:supremum(d.domain) for d in domains]
analytic_sol_func(t,x) =  sum([(8/(k^3*pi^3)) * sin(k*pi*x)*cos(C*k*pi*t) for k in 1:2:50000])

u_predict = reshape([first(phi([t,x],res.minimizer)) for t in ts for x in xs],(length(ts),length(xs)))
u_real = reshape([analytic_sol_func(t,x) for t in ts for x in xs], (length(ts),length(xs)))

diff_u = abs.(u_predict .- u_real)
p1 = plot(ts, xs, u_real, linetype=:contourf,title = "analytic");
p2 =plot(ts, xs, u_predict, linetype=:contourf,title = "predict");
p3 = plot(ts, xs, diff_u,linetype=:contourf,title = "error");
plot(p1,p2,p3)


text:
using JUGRNModelGenerator

# setup paths -
path_to_input_file = "./test/data/Model.json"
path_to_output_dir = "./test/generated_code"

# build the model object -
julia_model_object = build_julia_model_object(path_to_input_file, path_to_output_dir) |> check


text:
module FractionalDiffEq

using LinearAlgebra
using SpecialFunctions
using SparseArrays
using ApproxFun: Fun, Ultraspherical, coefficients
using InvertedIndices
using SpecialMatrices
using FFTW: fft, ifft
using UnPack: @unpack
using LoopVectorization: @turbo
using HypergeometricFunctions
using ToeplitzMatrices
using RecipesBase
using ForwardDiff

# Single-term fractional ordinary differential equations
include("singletermfode/PECE.jl")
include("singletermfode/PI.jl")
include("singletermfode/GL.jl")
include("singletermfode/ChebSpectral.jl")

# Multi-terms fractional ordinary differential equations
include("multitermsfode/matrix.jl")
include("multitermsfode/hankelmatrix.jl")
include("multitermsfode/ClosedForm.jl")
include("multitermsfode/highprecision.jl")
include("multitermsfode/MTPIEX.jl")
include("multitermsfode/MTPIPECE.jl")
include("multitermsfode/MTPIIMTrap.jl")
include("multitermsfode/MTPIIMRect.jl")

# System of fractional ordinary differential equations
include("fodesystem/PIPECE.jl")
include("fodesystem/FLMMBDF.jl")
include("fodesystem/FLMMNewtonGregory.jl")
include("fodesystem/FLMMTrap.jl")
include("fodesystem/PIEX.jl")
include("fodesystem/GLWithMemory.jl")
include("fodesystem/NonLinear.jl")

# Fractional partial differential equations
include("FPDE/FiniteDiffEx.jl")
include("FPDE/FiniteDiffIm.jl")
include("FPDE/CaputoDiscrete.jl")
include("FPDE/GL.jl")

# Fractional delay differential equations
include("FDDE/DelayPECE.jl")
include("FDDE/DelayPECEMultipleLags.jl")
include("FDDE/DelayPI.jl")
include("FDDE/Matrix.jl")
include("FDDE/DelayABM.jl")
include("FDDE/DelayABMSystem.jl")

# Distributed order differential equations
include("dode/utils.jl")
include("dode/matrix.jl")

# Fractional Differences equations
include("FractionalDifferences/PECE.jl")
include("FractionalDifferences/GL.jl")

# Fractional integral equations
include("FIE/Qmat.jl")
include("FIE/mycoeffs.jl")
include("FIE/main.jl")

# Mittag Leffler function
include("mlfun.jl")

include("utils.jl")
include("auxiliary.jl")

export FractionalDiffEqAlgorithm

# General types
export solve, FDEProblem

# FPDE problems
export FPDEProblem

# FDDE problems
export FDDEProblem, FDDESystem, FDDEMatrixProblem

# FODE problems
export SingleTermFODEProblem, MultiTermsFODEProblem, FODESystem, DODEProblem

# Fractional Difference probelms
export FractionalDifferenceProblem, FractionalDifferenceSystem

# FIE problems
export FIEProblem

###################################################

export AbstractFDESolution
export FODESolution, FIESolution, FDifferenceSolution, DODESolution

# FODE solvers
export PIEX, PIPECE, PIIMRect, PIIMTrap
export PECE, FODEMatrixDiscrete, ClosedForm, ClosedFormHankelM, ClosedFormHighPercision, GL
export ChebSpectral

# FPDE solvers
export FPDEMatrixDiscrete, FiniteDiffEx, FiniteDiffIm, ADV_DIF, GLDiff

# System of FODE solvers
export NonLinearAlg, GLWithMemory, FLMMBDF, FLMMNewtonGregory, FLMMTrap, PIEX

# FDDE solvers
export DelayPECE, DelayPI, DelayABM, MatrixForm

# DODE solvers
export DOMatrixDiscrete

# Fractional Differences Equations solvers
export PECEDifference

# Fractional Integral Equations solvers
export SpectralUltraspherical

###################################################

# Export some api to construct the equation
export eliminator, RieszMatrix, omega, meshgrid

# Export some special models
export bagleytorvik, diffusion

# Auxiliary functions
export mittleff, mittleffderiv

# Distributed order auxiliary SpecialFunctions
export DOB, DOF, DORANORT, isFunction

export myeval, ourfft, ourifft, rowfft, FastConv

end


text:
using LinearAlgebra, StatisticalNonlocality
import StatisticalNonlocality: ou_transition_matrix
import StatisticalNonlocality: fourier

filename = "two_state.jld2"
n = 1
M = ou_transition_matrix(n)

Λ, V = eigen(M)
Λ, W = eigen(M')

N = 128

D, x = fourier(N, a = 0, b = 2π)

γ = 1.0
U = 1.0
κ = 1.0

fluxkernel = -U * inv(κ .* D^2 - γ * I) * U
λ = eigvals(fluxkernel)
σ = svdvals(fluxkernel)

data_directory = "data"
mkpath(data_directory)
filepath = data_directory * "/" * filename
if isfile(filepath)
    rm(filepath)
end

file = jldopen(filepath, "a+")

file["kernel"] = fluxkernel
file["eigenvalues"] = λ
file["singularvalues"] = σ

close(file)

max_stars_repo_name: sandreza/StatisticalNonlocality

max_issues_repo_name: sandreza/StatisticalNonlocality

max_forks_repo_name: sandreza/StatisticalNonlocality


text:
using .UltraNest

export UltraReact

struct UltraReact <: ComradeSoss.AbstractNested end

"""
    dynesty_sampler(logj; nlive=400, kwargs...)
Takes a problem defined by the logj likelihood object and runs the
dynesty on it using the default options and static sampler.

Returns a chain, state
"""
function StatsBase.sample(::UltraReact, lj::Soss.ConditionalModel; kwargs...)
    lklhd, prt, tc, unflatten = ComradeSoss._split_conditional(lj)

    lklhdvec(X) = lklhd.(eachrow(X))
    vpnames = String["p$i" for i in 1:dimension(tc)]
    prtransvec(X) = reduce(vcat, (x -> prt(x)').(eachrow(X)))

    sampler = ultranest.ReactiveNestedSampler(
                                        vpnames,
                                        lklhdvec,
                                        transform=prtransvec,
                                        vectorized=true
                                        )
    res = sampler.run(;kwargs...)

    #res = sampler.results
    samples, weights = res["weighted_samples"]["points"], res["weighted_samples"]["weights"]
    logz = res["logz"][end]
    logzerr = res["logzerr"][end]
    logl = res["logl"][end]

    stats = (logz=logz, logzerr=logzerr, logl=logl)
    return ComradeSoss._create_tv(unflatten, samples, weights), stats
end


text:
"""
HelloDB saves data into files. Data is broken up into blocks, and
fetched/manipulated within a `Page` buffer.

Strings are saved withing a page as an int, representing the length of the
string, followed by chars.
"""
module FileMgrs
  using HelloDB
  using HelloDB.Blocks
  using HelloDB.Data: BLOCK_SIZE
  ##############################################################################
  # Exports
  ##############################################################################

  export FileMgr
  export append

  ##############################################################################
  # Implementation
  ##############################################################################

  type FileMgr
    dbdir::AbstractString
    isNew::Bool
    openfiles::Dict{AbstractString, IOStream}

    function FileMgr(dbname::AbstractString)
      # create db directory
      dbdir = homedir() * "/$dbname"
      isnew = false

      if !isdir(dbdir)
        isnew = true
        log("creating new dbdir file at=($(dbdir))")
        mkdir(dbdir)
        _make_dbdir_file(dbdir)
      end

      # clear leftover temp files
      tempfiles = filter(f->startswith(f, "temp"), readdir(dbdir))
      log("clearing tempfiles")
      [rm(dbdir + "/$file") for file in tempfiles]

      log("initializing new FileMgr")
      new(dbdir, isnew, Dict())
    end

    function FileMgr()
      new("foo", false, Dict())
    end
  end

  """
  Fill `buf` based on `block` filename and position.
  """
  function Base.read(filemgr::FileMgr, block::Block, buf::IOBuffer)
    _clear(buf)
    file = _getfile(filemgr, block.filename)
    seek(file, block.number * buf.maxsize)
    write(buf, readbytes(file, buf.maxsize))
  end

  function Base.write(filemgr::FileMgr, block::Block, buf::IOBuffer)
    seekstart(buf)
    file = _getfile(filemgr, block.filename)
    seek(file, block.number * buf.maxsize)
    written = write(file, readbytes(buf))
    flush(file)
    written
  end

  """
  Append contents of `buffer` to file `filename`.
  """
  function append(filemgr::FileMgr, filename::AbstractString, buffer::IOBuffer)
    file = _getfile(filemgr, filename)
    blocknum = size(filemgr, filename)
    block = Block(filename, blocknum)
    write(filemgr, block, buffer)
    block
  end

  """
  Size of file `filename`, in blocks`::Block`.
  """
  function size(filemgr::FileMgr, filename::AbstractString)
    file = _getfile(filemgr, filename)
    ratio = filesize(file) / BLOCK_SIZE
    if ratio % BLOCK_SIZE != 0
      round(Int, floor(ratio)) + 1
    else
      ratio
    end
  end

  ##############################################################################
  # Private
  ##############################################################################

  function _clear(buf::IOBuffer)
    buf.size = 0
  end

  # TODO: better filename handling.
  function _getfile(filemgr::FileMgr, filename)
    filename = basename(filename)
    if haskey(filemgr.openfiles, filename)
      filemgr.openfiles[filename]
    else
      file = open(filemgr.dbdir * "/$filename", "a+")
      filemgr.openfiles[filename] = file
      file
    end
  end

  """
  Will prevent you from catastrophe.
  """
  function _make_dbdir_file(dbdir)
    filepath = dbdir * "/.dbdir"
    file = open(dbdir * "/.dbdir", "w")
    close(file)
  end
end


text:
using Documenter, Microbiome, Microbiome.Dictionaries

makedocs(
    sitename = "Microbiome.jl",
    strict=:doctest,
    pages = [
        "Home" => "index.md",
        "Samples and features" => "samples_features.md",
        "Profiles and Communities" => "profiles.md",
        "Diversity measures" => "diversity.md",
    ],
    authors = "Kevin Bonham, PhD <kbonham@wellesley.edu",
    format = Documenter.HTML(
        prettyurls = get(ENV, "CI", nothing) == "true",
        edit_link="main",
        canonical="http://docs.ecojulia.org/BiobakeryUtils.jl/stable/")
        )

deploydocs(
    repo = "github.com/EcoJulia/Microbiome.jl.git",
    push_preview=true,
    devbranch="main"
)


text:
function _clean_file(filename::String)
    coverage_filename = string(filename, ".microcov")
    rm(coverage_filename; force = true, recursive = true)
    return nothing
end

function _clean_directory(path::String)
    microcov_file_suffix = ".microcov"
    for (root, dirs, files) in walkdir(".")
        for file in files
            if endswith(file, microcov_file_suffix)
                rm(joinpath(root, file); force = true, recursive = true)
            end
        end
    end
    return nothing
end


text:
function mark_around_border!(r)
    height = 11
    width = 12
    side1 = Sud
    side2 = West
    go_to_corner_with_walls!(r, side1, side2)

    k = 0
    number_now = 0
    side = Ost
    while !isborder(r, Nord)
        if k % 2 == 0
            number_now = counter(Ost)
            side = Ost
        else
            number_now = counter(West)
            side = West
        end
        k += 1
        if (number_now + 1) != width
            break
        end
        move!(r, Nord)
    end

    if (side == West)
        for i in [Sud, West, Nord, Ost, Sud]
            mark_around(i)
        end
    else
        for i in [Nord, Ost, Sud, West, Nord]
            mark_around(i)
        end
    end

    while !isborder(r, side1) || !isborder(r, side2)
        while (!isborder(r, side1))
            move!(r, side1)
        end
        while (!isborder(r, side2))
            move!(r, side2)
        end
    end

    while length(sides) > 0
        move!(r, inverse_side(pop!(sides)))
    end
end

max_stars_repo_name: dmitrii-a-piskunov/mirea-progs

max_issues_repo_name: dmitrii-a-piskunov/mirea-progs

max_forks_repo_name: dmitrii-a-piskunov/mirea-progs