text:
"""
Package: Forecast

    function d(x::{AbstractVector,AbstractArray,TimeArray},
               order::Int=1;
               lag::Int=1,
               center::Bool=true,
               pad::Bool=true)

Return Lagged differences of a given order for Vector, Array and TimeSeries.

# Arguments
- `x`: Vector or Array of data.
- `order`: Order of the differences; number of recursive iterations on the same vector/array.
- `lag`: Lag for the difference.
- `center`: Center the result in the response using Missing values.
- `pad`: Includes or removes `missing` pad.

# Returns
Lagged differences Vector or Array of a given order.

# Examples
```julia-repl
julia> x = [1,2,3,4,5];
julia> d(x)
5-element Array{Any,1}:
 1
 1
 1
 1
  missing

julia> d(x,2)
5-element Array{Any,1}:
  missing
 0
 0
 0
  missing

julia> d(x;lag=2,pad=false)
3-element Array{Any,1}:
 2
 2
 2

julia> x = reshape(collect(1:20),10,2)
10×2 Array{Int64,2}:
  1  11
  2  12
  3  13
  4  14
  5  15
  6  16
  7  17
  8  18
  9  19
 10  20

julia> d(x,2;lag=2,pad=false)
6×2 Array{Any,2}:
 0  0
 0  0
 0  0
 0  0
 0  0
 0  0

julia> d(co2())
[ Info: Dataset used in Cleveland et al. paper
4609×1 TimeArray{Any,1,Date,Array{Any,1}} 1974-05-17 to 1986-12-31
│            │ A       │
├────────────┼─────────┤
│ 1974-05-17 │ -0.27   │
│ 1974-05-18 │ 0.35    │
│ 1974-05-19 │ 0.18    │
   ⋮
│ 1986-12-30 │ missing │
│ 1986-12-31 │ missing │

```
"""
function d(x::AbstractVector,
           order::Int=1;
           lag::Int=1,
           center::Bool=true,
           pad::Bool=true)

    N = length(x)
    @assert 0 <= lag & lag <= (N-1)

    if (lag == 0) | (order == 0)
        return x
    end

    dx = Vector{Any}(missing,N)

    # update dx function
    function udx(v)
        if (lag == 1)
            dx[1:(N-lag)] = diff(v)
        else
            dx[1:(N-lag)] = diff(hcat(circshift(v,lag),v),dims=2)[(1+lag):end]
        end
    end

    # update dx with x
    udx(x)

    # update dx with previous dx
    for i in 2:order
        udx(dx)
    end

    # center values
    if pad
        return circshift(dx,div(lag,2)+div(order,2))
    else
        ivp = findfirst(!ismissing, dx)
        fvp = findlast(!ismissing, dx)
        return(dx[ivp:fvp])
    end

end

function d(x::AbstractArray,
           order::Int=1;
           lag::Int=1,
           center::Bool=true,
           pad::Bool=true)

    N,M = size(x)

    @assert 0 <= lag & lag <= (N-1)

    if (lag == 0) | (order == 0)
        return x
    end

    dx = Array{Any}(missing,N,M)

    # update dx function
    function udx(v)
        if (lag == 1)
            dx[1:(N-lag),:] = diff(v,dims=1)
        else
            xh = hcat(circshift(v,lag),v)
            dx[1:(N-lag),:] = xh[(1+lag):end,M+1:2*M] - xh[(1+lag):end,1:M]
        end
    end

    # update dx with x
    udx(x)

    # update dx with previous dx
    for i in 2:order
        udx(dx)
    end

    # center values
    if pad
        return circshift(dx,div(lag,2)+div(order,2))
    else
        ivp = findfirst(!ismissing, dx)
        fvp = findlast(!ismissing, dx)
        return(dx[ivp:fvp])
    end

end

function d(x::TimeArray,
           order::Int=1;
           lag::Int=1,
           center::Bool=true,
           pad::Bool=true)

    vx = values(x)
    tsx = timestamp(x)
    replace!(vx, NaN => missing)

    dvx = d(vx, order; lag=lag, center=center, pad=pad)

    TimeArray(tsx,dvx)

end


text:
"""
    tex_layout(mathexpr::TeXExpr, font_family)

Recursively determine the layout of the math expression represented the given
TeXExpr for the given font set.
"""
tex_layout(expr, font_family::FontFamily) = tex_layout(expr, LayoutState(font_family))

function tex_layout(expr, state)
    font_family = state.font_family
    head = expr.head
    args = [expr.args...]
    shrink = 0.6

    try
        if head in [:char, :delimiter, :digit, :punctuation, :symbol]
            char = args[1]
            return TeXChar(char, state, head)
        elseif head == :combining_accent
            accent, core = tex_layout.(args, state)

            y = topinkbound(core) - xheight(font_family)

            if core.slanted
                α = slant_angle(font_family)
                x = (y + bottominkbound(accent)) * tan(α) / 2
            else
                x = 0.0
            end

            return Group(
                [core, accent],
                Point2f[
                    (0, 0),
                    (x + hmid(core) - hmid(accent), y)
                ],
                [1, 1]
            )
        elseif head == :decorated
            core, sub, super = tex_layout.(args, state)

            core_width = advance(core)

            return Group(
                [core, sub, super],
                Point2f[
                    (0, 0),
                    (core_width, -0.2),
                    (core_width, xheight(core) - 0.5 * descender(super))],
                [1, shrink, shrink]
            )
        elseif head == :delimited
            elements = tex_layout.(args, state)
            left, content, right = elements

            height = inkheight(content)
            left_scale = max(1, height / inkheight(left))
            right_scale = max(1, height / inkheight(right))
            scales = [left_scale, 1, right_scale]

            dxs = advance.(elements) .* scales
            xs = [0, cumsum(dxs[1:end-1])...]

            # TODO Height calculation for the parenthesis looks wrong
            # TODO Check what the algorithm should be there
            # Center the delimiters in the middle of the bot and top baselines ?
            return Group(elements,
                Point2f[
                    (xs[1], -bottominkbound(left) + bottominkbound(content)),
                    (xs[2], 0),
                    (xs[3], -bottominkbound(right) + bottominkbound(content))
                ],
                scales
            )
        elseif head == :font
            modifier, content = args
            return tex_layout(content, add_font_modifier(state, modifier))
        elseif head == :frac
            numerator = tex_layout(args[1], state)
            denominator = tex_layout(args[2], state)

            # extend fraction line by half an xheight
            xh = xheight(font_family)
            w = max(inkwidth(numerator), inkwidth(denominator)) + xh/2

            # fixed width fraction line
            lw = thickness(font_family)

            line = HLine(w, lw)
            y0 = xh/2 - lw/2

            # horizontal center align for numerator and denominator
            x1 = (w-inkwidth(numerator))/2
            x2 = (w-inkwidth(denominator))/2

            ytop    = y0 + xh/2 - bottominkbound(numerator)
            ybottom = y0 - xh/2 - topinkbound(denominator)

            return Group(
                [line, numerator, denominator],
                Point2f[(0,y0), (x1, ytop), (x2, ybottom)]
            )
        elseif head == :function
            name = args[1]
            elements = TeXChar.(collect(name), state, Ref(:function))
            return horizontal_layout(elements)
        elseif head == :group || head == :expr
            elements = tex_layout.(args, state)
            return horizontal_layout(elements)
        elseif head == :integral
            pad = 0.1
            sub, super = tex_layout.(args[2:3], state)

            # Always use ComputerModern fallback for the integral sign
            # as the Unicode LaTeX approach requires to use glyph variant
            # which is unlikely to be supported by backends
            intfont = load_font(joinpath("ComputerModern", "cmex10.ttf"))
            int = TeXChar(Char(0x5a), intfont)
            h = inkheight(int)

            return Group(
                [int, sub, super],
                Point2f[
                    (0, h/2 + xheight(font_family)/2),
                    (
                        0.15 - inkwidth(sub)*shrink/2,
                        -h/2 + xheight(font_family)/2 - topinkbound(sub)*shrink - pad
                    ),
                    (
                        0.85 - inkwidth(super)*shrink/2,
                        h/2 + xheight(font_family)/2 + pad
                    )
                ],
                [1, shrink, shrink]
            )
        elseif head == :space
            return Space(args[1])
        elseif head == :spaced
            sym = tex_layout(args[1], state)
            return horizontal_layout([Space(0.2), sym, Space(0.2)])
        elseif head == :sqrt
            content = tex_layout(args[1], state)
            sqrt = TeXChar('√', state, :symbol)

            relpad = 0.15

            h = inkheight(content)
            ypad = relpad * h
            h += 2ypad

            if h > inkheight(sqrt)
                sqrt = TeXChar('⎷', state, :symbol)
            end

            h = max(inkheight(sqrt), h)

            # The root symbol must be manually placed
            y0 = bottominkbound(content) - bottominkbound(sqrt) - ypad/2
            y = y0 + bottominkbound(sqrt) + h
            xpad = advance(sqrt) - inkwidth(sqrt)
            w =  inkwidth(content) + 2xpad

            lw = thickness(font_family)
            hline = HLine(w, lw)
            vline = VLine(inkheight(sqrt) - h, lw)

            return Group(
                [sqrt, vline, hline, content],
                Point2f[
                    (0, y0),
                    (rightinkbound(sqrt) - lw/2, y),
                    (rightinkbound(sqrt) - lw/2, y - lw/2),
                    (advance(sqrt), 0)
                ]
            )

        elseif head == :underover
            core, sub, super = tex_layout.(args, state)

            mid = hmid(core)
            dxsub = mid - hmid(sub) * shrink
            dxsuper = mid - hmid(super) * shrink

            under_offset = bottominkbound(core) - (ascender(sub) - xheight(sub)/2) * shrink
            over_offset = topinkbound(core) - descender(super)

            # The leftmost element must have x = 0
            x0 = -min(0, dxsub, dxsuper)

            return Group(
                [core, sub, super],
                Point2f[
                    (x0, 0),
                    (x0 + dxsub, under_offset),
                    (x0 + dxsuper, over_offset)
                ],
                [1, shrink, shrink]
            )
        end
    catch
        # TODO Better error
        rethrow()
        @error "Error while processing expr"
    end

    @error "Unsupported head $(head) in expr:\n$expr"
end

tex_layout(::Nothing, state) = Space(0)

"""
    horizontal_layout(elements)

Layout the elements horizontally, like normal text.
"""
function horizontal_layout(elements)
    dxs = advance.(elements)
    xs = [0, cumsum(dxs[1:end-1])...]

    return Group(elements, Point2f.(xs, 0))
end

function layout_text(string, font_family)
    isempty(string) && return Space(0)

    elements = TeXChar.(collect(string), LayoutState(font_family), Ref(:text))
    return horizontal_layout(elements)
end

"""
    unravel(element::TeXElement, pos, scale)

Flatten the layouted TeXElement and produce a single list of base element with
their associated absolute position and scale.
"""
function unravel(group::Group, parent_pos=Point2f(0), parent_scale=1.0f0)
    scales = group.scales .* parent_scale
    positions = [parent_pos .+ pos for pos in parent_scale .* group.positions]
    elements = []

    for (elem, pos, scale) in zip(group.elements, positions, scales)
        push!(elements, unravel(elem, pos, scale)...)
    end

    return elements
end

unravel(::Space, pos, scale) = []
unravel(element, pos, scale) = [(element, pos, scale)]

"""
    generate_tex_elements(str)

Create a list of tuple `(texelement, position, scale)` from a string
of LaTeX math mode code. The elements' positions and scales are such as to
approximatively reproduce the LaTeX output.

The elments are of one of the following types
    - `TeXChar` a (unicode) character, in a specific font.
    - `HLine` a horizontal line.
    - `VLine` a vertical line.
"""
function generate_tex_elements(str, font_family=FontFamily())
    expr = texparse(str)
    layout = tex_layout(expr, font_family)
    return unravel(layout)
end

# Still hacky as hell
function generate_tex_elements(str::LaTeXString, font_family=FontFamily())
    parts = String.(split(str, raw"$"))
    groups = Vector{TeXElement}(undef, length(parts))
    texts = parts[1:2:end]
    maths = parts[2:2:end]

    groups[1:2:end] = layout_text.(texts, Ref(font_family))
    groups[2:2:end] = tex_layout.(texparse.(maths), Ref(font_family))

    return unravel(horizontal_layout(groups))
end

ext: jl


text:
using Distributions, Plots, LaTeXStrings; pyplot()

function statPair(dist,n)
    sample = rand(dist,n)
    [mean(sample),var(sample)]
end

stdUni = Uniform(-sqrt(3),sqrt(3))
n, N = 3, 10^5

dataUni     = [statPair(stdUni,n) for _ in 1:N]
dataUniInd  = [[mean(rand(stdUni,n)),var(rand(stdUni,n))] for _ in 1:N]
dataNorm    = [statPair(Normal(),n) for _ in 1:N]
dataNormInd = [[mean(rand(Normal(),n)),var(rand(Normal(),n))] for _ in 1:N]

p1 = scatter(first.(dataUni), last.(dataUni),
	c=:blue, ms=1, msw=0, label="Same group")
p1 = scatter!(first.(dataUniInd), last.(dataUniInd),
	c=:red, ms=0.8, msw=0, label="Separate group",
	xlabel=L"\overline{X}", ylabel=L"S^2")

p2 = scatter(first.(dataNorm), last.(dataNorm),
	c=:blue, ms=1, msw=0, label="Same group")
p2 = scatter!(first.(dataNormInd), last.(dataNormInd),
	c=:red, ms=0.8, msw=0, label="Separate group",
	xlabel=L"\overline{X}", ylabel=L"$S^2$")

plot(p1, p2, ylims=(0,5), size=(800, 400))

max_stars_repo_name: Yoshinobu-Ishizaki/StatsWithJuliaBook

max_issues_repo_name: Yoshinobu-Ishizaki/StatsWithJuliaBook

max_forks_repo_name: Yoshinobu-Ishizaki/StatsWithJuliaBook


text:
using NumericalMethodsforEngineers
using Test

af = Float64[16, 4, 5, 8, -4, 22]
ai = [16, 4, 5, 8, -4, 22]
ac = [16+2im, 4-2im, 5+0im, 8, -4-4im, 22+9im]
acf = [16.0+2im, 4-2.0im, 5+0.0im, 8.0, -4-4im, 22.0+9im]

kdiag =[1, 3, 6]

sres=SkylineMatrix(kdiag, af)

@test full(sres) == [
  16.0   4.0   8.0;
   4.0   5.0  -4.0;
   8.0  -4.0  22.0
]

@test full(sres) == fromskyline(kdiag, af)

sres=SkylineMatrix(kdiag, ai)
@test full(sres) == fromskyline(kdiag, ai)

sres=SkylineMatrix(kdiag, ac)
@test full(sres) == fromskyline(kdiag, ac)

sres=SkylineMatrix(kdiag, acf)
@test full(sres) == fromskyline(kdiag, acf)


text:
# Operating at these extreme ocean depths has overloaded the submarine's
# reactor; it needs to be rebooted.
#
# The reactor core is made up of a large 3-dimensional grid made up entirely of
# cubes, one cube per integer 3-dimensional coordinate (x,y,z). Each cube can
# be either on or off; at the start of the reboot process, they are all off.
# (Could it be an old model of a reactor you've seen before?)
#
# To reboot the reactor, you just need to set all of the cubes to either on or
# off by following a list of reboot steps (your puzzle input). Each step
# specifies a cuboid (the set of all cubes that have coordinates which fall
# within ranges for x, y, and z) and whether to turn all of the cubes in that
# cuboid on or off.
#
# For example, given these reboot steps:
#
# on x=10..12,y=10..12,z=10..12
# on x=11..13,y=11..13,z=11..13
# off x=9..11,y=9..11,z=9..11
# on x=10..10,y=10..10,z=10..10
#
# The first step (on x=10..12,y=10..12,z=10..12) turns on a 3x3x3 cuboid
# consisting of 27 cubes:
#
# 10,10,10
# 10,10,11
# 10,10,12
# 10,11,10
# 10,11,11
# 10,11,12
# 10,12,10
# 10,12,11
# 10,12,12
# 11,10,10
# 11,10,11
# 11,10,12
# 11,11,10
# 11,11,11
# 11,11,12
# 11,12,10
# 11,12,11
# 11,12,12
# 12,10,10
# 12,10,11
# 12,10,12
# 12,11,10
# 12,11,11
# 12,11,12
# 12,12,10
# 12,12,11
# 12,12,12
#
# The second step (on x=11..13,y=11..13,z=11..13) turns on a 3x3x3 cuboid that
# overlaps with the first. As a result, only 19 additional cubes turn on; the
# rest are already on from the previous step:
#
# 11,11,13
# 11,12,13
# 11,13,11
# 11,13,12
# 11,13,13
# 12,11,13
# 12,12,13
# 12,13,11
# 12,13,12
# 12,13,13
# 13,11,11
# 13,11,12
# 13,11,13
# 13,12,11
# 13,12,12
# 13,12,13
# 13,13,11
# 13,13,12
# 13,13,13
#
# The third step (off x=9..11,y=9..11,z=9..11) turns off a 3x3x3 cuboid that
# overlaps partially with some cubes that are on, ultimately turning off 8
# cubes:
#
# 10,10,10
# 10,10,11
# 10,11,10
# 10,11,11
# 11,10,10
# 11,10,11
# 11,11,10
# 11,11,11
#
# The final step (on x=10..10,y=10..10,z=10..10) turns on a single cube,
# 10,10,10. After this last step, 39 cubes are on.
#
# The initialization procedure only uses cubes that have x, y, and z positions
# of at least -50 and at most 50. For now, ignore cubes outside this region.
#
# Here is a larger example:
#
# on x=-20..26,y=-36..17,z=-47..7
# on x=-20..33,y=-21..23,z=-26..28
# on x=-22..28,y=-29..23,z=-38..16
# on x=-46..7,y=-6..46,z=-50..-1
# on x=-49..1,y=-3..46,z=-24..28
# on x=2..47,y=-22..22,z=-23..27
# on x=-27..23,y=-28..26,z=-21..29
# on x=-39..5,y=-6..47,z=-3..44
# on x=-30..21,y=-8..43,z=-13..34
# on x=-22..26,y=-27..20,z=-29..19
# off x=-48..-32,y=26..41,z=-47..-37
# on x=-12..35,y=6..50,z=-50..-2
# off x=-48..-32,y=-32..-16,z=-15..-5
# on x=-18..26,y=-33..15,z=-7..46
# off x=-40..-22,y=-38..-28,z=23..41
# on x=-16..35,y=-41..10,z=-47..6
# off x=-32..-23,y=11..30,z=-14..3
# on x=-49..-5,y=-3..45,z=-29..18
# off x=18..30,y=-20..-8,z=-3..13
# on x=-41..9,y=-7..43,z=-33..15
# on x=-54112..-39298,y=-85059..-49293,z=-27449..7877
# on x=967..23432,y=45373..81175,z=27513..53682
#
# The last two steps are fully outside the initialization procedure area; all
# other steps are fully within it. After executing these steps in the
# initialization procedure region, 590784 cubes are on.
#
# Execute the reboot steps. Afterward, considering only cubes in the region
# x=-50..50,y=-50..50,z=-50..50, how many cubes are on?

const (minx,maxx) = (-50,50)
const (miny,maxy) = (-50,50)
const (minz,maxz) = (-50,50)

# build grid - read input
grid = falses(maxz-minz+1, maxy-miny+1, maxx-minx+1)
foreach(eachline("input.txt")) do l
  # parse line
  (instr,coords...) = match(
    r"(on|off) x=(-?\d+)..(-?\d+),y=(-?\d+)..(-?\d+),z=(-?\d+)..(-?\d+)",
    l
  )
  (x1,x2,y1,y2,z1,z2) = parse.(Int, coords)

  # ensure the ranges cross inside the initialization area
  all([
    minx <= x1 <= maxx || minx <= x2 <= maxx,
    miny <= y1 <= maxy || miny <= y2 <= maxy,
    minz <= z1 <= maxz || minz <= z2 <= maxz
  ]) || return

  # clamp and translate
  x1 = clamp(x1, minx, maxx) - minx + 1
  x2 = clamp(x2, minx, maxx) - minx + 1
  y1 = clamp(y1, miny, maxy) - miny + 1
  y2 = clamp(y2, miny, maxy) - miny + 1
  z1 = clamp(z1, minz, maxz) - minz + 1
  z2 = clamp(z2, minz, maxz) - minz + 1

  # update
  grid[z1:z2, y1:y2, x1:x2] .= instr == "on"
end

println("Result: ", count(grid))


text:
module DualMatrixTools

using DualNumbers, LinearAlgebra, SparseArrays, SuiteSparse

# Personal add-on: Allow for multiplication by ε of a sparse matrix
import Base.*
"""
    *(d::Dual, M::SparseMatrixCSC)

Exactly multiplies a sparse matrix `M` by the dual number `d`.
This overload allows `ε * M` to be non-zero.
(This is to correct the effect of the `ε == 0` feature in DualNumbers.jl,
which silently removes non-zeros from sparse matrices when multiplied by `ε`.)
"""
function *(d::Dual, M::SparseMatrixCSC)
    i, j, v = findnz(M)
    m, n = size(M)
    return sparse(i, j, d .* v, m, n)
end
export *

"""
    DualFactors

Container type to work efficiently with backslash on dual-valued sparse matrices.

`factorize(M)` will create an instance containing
- `Af = factorize(realpart.(M))` — the factors of the real part
- `B = dualpart.(M)` — the dual part
for a dual-valued matrix `M`.

This is because only the factors of the real part are needed when solving a linear system of the type ``M x = b`` for a dual-valued matrix ``M = A + \\varepsilon B``.
In fact, the inverse of ``M`` is given by
``M^{-1} = (I - \\varepsilon A^{-1} B) A^{-1}``.
"""
mutable struct DualFactors{TAf,TB}
    Af::TAf # the factors of the real part
    B ::TB  # the ε part
end
export DualFactors

# Factorization functions
for f in (:lu, :qr, :cholesky, :factorize)
    @eval begin
        import LinearAlgebra: $f
        """
        $($f)(M::SparseMatrixCSC{<:Dual,<:Int})

        Invokes `$($f)` on just the real part of `M` and stores it along with the dual part into a `DualFactors` object.
        """
        $f(M::SparseMatrixCSC{<:Dual,<:Int}) = DualFactors($f(realpart.(M)), dualpart.(M))
        """
        $($f)(M::Array{<:Dual,2})

        Invokes `$($f)` on just the real part of `M` and stores it along with the dual part into a `DualFactors` object.
        """
        $f(M::Array{<:Dual,2}) = DualFactors($f(realpart.(M)), dualpart.(M))
        export $f
    end
end

# In-place factorization for the case where the real part is already stored
import LinearAlgebra: factorize
function factorize(Mf::DualFactors, M; update_factors = false)
    Mf.B = dualpart.(M)
    if update_factors
        Mf.Af = factorize(realpart.(M))
    end
    return Mf
end
export factorize

# Adjoint and transpose definitions for `DualFactors`
for f in (:adjoint, :transpose)
    @eval begin
        import Base: $f
        """
        $($f)(M::DualFactors)

        Invokes `$($f)` on both `M.Af` and `M.B` and returns them into a new `DualFactors` object.
        """
        $f(M::DualFactors) = DualFactors($f(M.Af), $f(M.B))
        export $f
    end
end

import Base.\
"""
    \\(M::DualFactors, y::AbstractVecOrMat{Float64})

Backsubstitution for `DualFactors`.
See `DualFactors` for details.
"""
function \(M::DualFactors, y::AbstractVecOrMat{Float64})
    A, B = M.Af, M.B
    A⁻¹y = A \ y
    return A⁻¹y - ε * (A \ (B * A⁻¹y))
end

"""
    \\(M::DualFactors, y::AbstractVecOrMat{Dual128})

Backsubstitution for `DualFactors`.
See `DualFactors` for details.
"""
function \(M::DualFactors, y::AbstractVecOrMat{Dual128})
    a, b = realpart.(y), dualpart.(y)
    A, B = M.Af, M.B
    A⁻¹a = A \ a
    return A⁻¹a + ε * (A \ (b - B * A⁻¹a))
end

"""
    \\(Af::Factorization{Float64}, y::AbstractVecOrMat{Dual128})

Backsubstitution for Dual-valued RHS.
"""
function \(Af::Factorization{Float64}, y::AbstractVecOrMat{Dual128})
    return (Af \ realpart.(y)) + ε * (Af \ dualpart.(y))
end

"""
    \\(M::AbstractArray{Dual128,2}, y::AbstractVecOrMat)

Backslash (factorization and backsubstitution) for Dual-valued matrix `M`.
"""
\(M::SparseMatrixCSC{<:Dual,<:Int}, y::AbstractVecOrMat) = factorize(M) \ y
\(M::Array{<:Dual,2}, y::AbstractVecOrMat) = factorize(M) \ y
export \

import Base.isapprox
function isapprox(x::AbstractVecOrMat{Dual128}, y::AbstractVecOrMat{Dual128})
    bigx = [realpart.(x) dualpart.(x)]
    bigy = [realpart.(y) dualpart.(y)]
    return isapprox(bigx, bigy)
end
isapprox(x::AbstractVecOrMat, y::AbstractVecOrMat{Dual128}) = isapprox(dual.(x), y)
isapprox(x::AbstractVecOrMat{Dual128}, y::AbstractVecOrMat) = isapprox(x, dual.(y))
function isapprox(x::Dual128, y::Dual128)
    bigx = [realpart(x) dualpart(x)]
    bigy = [realpart(y) dualpart(y)]
    return isapprox(bigx, bigy)
end
isapprox(x::Float64, y::Dual128) = isapprox(dual(x), y)
isapprox(x::Dual128, y::Float64) = isapprox(x, dual(y))

end # module


text:
include("initialization_step.jl")
include("partition_assignment_step.jl")
include("codebook_update_step.jl")

"""
    function quantizer(data::AbstractMatrix, η::Loss, codebook::Any, config::NamedTuple)

Main Quantization function calling the Anisotorpic or L2 Product/Vector quantizer.
Takes an input data set, a loss function (L2 or Anisotropic Weights), an input codebook to start from
in case a starting estimate is provided and configuration settings. Returns a QuantizerData data structure.
"""
function quantizer(data::AbstractMatrix, η::L2_loss, codebook::Any, config::NamedTuple)
    if config[:verbose] @info("Starting L2 Quantization...") end
    thread, update_method, optim_method, processing = generate_config_vars(config)
    qd = GenerateQuantizerData(data, config.n_codebooks, config.n_centers)
    initialization!(data, qd, thread, codebook)
    dist = 1
    i = 0
    while (dist > config[:stopcond]) && (i < config[:max_iter])
        i += 1
        C_old = deepcopy(qd.C)
        iteration_loss = assignment_step!(data, qd, update_method, thread)
        update_codebook!(data, qd, update_method, processing, optim_method)
        dist = sqeuclidean(C_old, qd.C)
        if config[:verbose] && config[:multithreading] @info("Update distance: $dist")
        elseif config[:verbose] @info("Iteration loss: $iteration_loss\tUpdate distance:$dist") end
    end
    if config[:verbose]
        if i < config[:max_iter] @info("Product quantization converged after $i iterations")
        else @info("Training stopped after reaching max ($i) iterations") end
    end
    if config[:inverted_index] rebuild_Bmatrix!(qd) end
    return qd
end

function quantizer(data::AbstractMatrix, η::AnisotropicWeights, codebook::Any, config::NamedTuple)
    if config[:verbose] @info("Starting Anisotropic Quantization...") end
    thread, update_method, optim_method, processing = generate_config_vars(config)
    qd = GenerateQuantizerData(data, config.n_codebooks, config.n_centers)
    initialization!(data, qd, thread, codebook)
    assignment_step!(data, qd, AHPQ.BMatrix(), thread)
    dist = 1
    i = 0
    while (dist > config[:stopcond]) && (i < config[:max_iter])
        i += 1
        C_old = deepcopy(qd.C)
        iteration_loss = assignment_step!(data, qd, η, config[:max_iter_assignments], thread)
        update_codebook!(data, qd, η, processing, optim_method)
        dist = sqeuclidean(C_old, qd.C)
        if config[:verbose] && config[:multithreading] @info("Update distance: $dist")
        elseif config[:verbose] @info("Iteration loss: $iteration_loss\tUpdate distance:$dist") end
    end
    if config[:verbose]
        if i < config[:max_iter] @info("Product quantization converged after $i iterations")
        else @info("Training stopped after reaching max ($i) iterations") end
    end
    return qd
end

function logrange(x1, n) return (Int(floor(10^y)) for y in range(log10(x1)-n+1, log10(x1), length=n)) end
function subsample(n_dp, data) return data[:,shuffle!(collect(1:size(data)[2]))[1:n_dp]] end

"""
    function incremental_quantization(data::AbstractMatrix, η::Weights, config::NamedTuple)

Wrapper for the main quantization function to instead train the VQ/PQ on an increasingly large subsample
of the data set. Can be used to speed up training of the quantizer in case of high-dimensional data.
"""
function incremental_quantization(data::AbstractMatrix, η::Weights, config::NamedTuple)
    thread, update_method, optim_method, processing = generate_config_vars(config)
    tp = if config.training_points > 0 config.training_points else size(data)[2] end
    samplesizes = collect(logrange(tp, config[:increment_steps]))
    if config[:verbose] @info("\n\nFitting quantizer on $(samplesizes[1]) data points...") end
    qd = quantizer(subsample(samplesizes[1], data), η, 0, config)
    codebook = deepcopy(qd.C)

    for samplesize in samplesizes[2:end]
        if config[:verbose] @info("\n\nFitting quantizer on $samplesize data points...") end
        qd = quantizer(subsample(samplesize, data), η, codebook, merge(config, (;use_precomputed_codebook=true)))
        codebook = deepcopy(qd.C)
    end
    return qd
end


text:
using DrWatson
@quickactivate "StatReth"

# %%
using StatsBase
using Distributions
using Plots

# %% 3.20 - 3.23
pdf.(Binomial(2, 0.7), 0:2)

rand(Binomial(2, 0.7), 1)

rand(Binomial(2, 0.7), 10)

N = 10_000
dummy_w = rand(Binomial(2, 0.7), N)
counts(dummy_w) ./ N |> display   # either
countmap(dummy_w)                 # or, depending on what you want

# %% 3.24, 3.25
dummy_w = rand(Binomial(9, 0.7), N)
histogram(dummy_w)

w = rand(Binomial(9, 0.6), N)

# %% 3.26
samples = let
    # %% 3.2 - 3.5
    n = 1000
    p_grid = range(0, 1, length = n)
    prob_p = ones(n)
    prob_data = @. pdf(Binomial(9, p_grid), 6)
    posterior = prob_data .* prob_p
    posterior ./= sum(posterior)
    weights = pweights(posterior)
    sample(p_grid, weights, 10_000)
end

# Version 1: Imagine this as having all the binomial distributions from the p of the
# samples added up. The random numbers are drawn from the resulting distribution.
mix = MixtureModel(Binomial.(9, samples))
w1 = rand(mix, 10_000)

# Version 2: From each p in the sample you make a binomial distribution and draw 1 value
# from it. Then you put them all into one big array.
w2 = rand.(Binomial.(9, samples), 1)
w2 = vcat(w2...)

# These are two different approaches but give the same result for lots of values drawn.
# However, if you are only drawing a few values (less than samples) the first version
# is preferrable since you can still use all the information in the samples.

histogram(w1, normalize = :probability, alpha = 0.6)
histogram!(w2, normalize = :probability, alpha = 0.6)


text:
using CUDA, LinearAlgebra,Logging

"""
https://stats.stackexchange.com/questions/503058/relationship-between-cholesky-decomposition-and-matrix-inversion
https://www.lume.ufrgs.br/bitstream/handle/10183/151001/001009773.pdf
experiments with cholesky and inverse

code taken from
https://www.codesansar.com/numerical-methods/python-program-inverse-matrix-using-gauss-jordan.htm

and

https://discourse.julialang.org/t/cholesky-simple-example/69480/7

most important mathemathically
https://stats.stackexchange.com/questions/147210/efficient-fast-mahalanobis-distance-computation/147222#147222

"""

using LinearAlgebra
x = [4 5 9]
A = [ 6 15 55;15 55 225;55 225 979  ]

x* inv(A)* transpose(x)

L = zeros(3,3)
lower =L
matrix = A
n=3

for i in 1:n
    for j in 1:i
        sum1 = 0;
        # summation for diagonals
        if (j == i)
            for k in 1:j
                sum1 += (lower[j,k])^2
            end#for
            lower[j,j] = sqrt(matrix[j,j] - sum1)
        else
            # Evaluating L(i, j)  using L(j, j)
            for k in 1:j
                sum1 += (lower[i,k] * lower[j,k]);
            end#for
            if(lower[j,j] > 0)
                lower[i,j] = (matrix[i,j] - sum1) /lower[j,j]
            end#if
        end#if else
    end#for
end#for

L
cholesky([ 6 15 55;15 55 225;55 225 979  ]).L
# so now according tohttps://stats.stackexchange.com/questions/147210/efficient-fast-mahalanobis-distance-computation/147222#147222
# we need to calculate Ly=x by forward substitution example(https://www.youtube.com/watch?v=8Z5oy80po-A) where y is unknown
# and mahalanobis will be just y times y transpose so we square each entry and add them up

#basic
y=[0.0,0.0,0.0]
x = [4.0,5.0,9.0]
y[1]= x[1]/L[1,1]
y[2] = (x[2]-L[2,1]*y[1])/L[2,2]
#y[3]= (x[3] -y[2]*L[3,2] - y[3]*L[3,3])/L[3,1]
#y[1]* L[3,1]+ y[2]*L[3,2] + y[3]*L[3,3]= x[3]
y[3]= (x[3]-y[2]*L[3,2]-y[1]* L[3,1])/L[3,3]

L*y
x

######## correct main equation
using LinearAlgebra
x = [4.0,5.0,9.0]
A = [ 6.0 15.0 55.0;
    15.0 55.0 225.0;
    55.0 225.0 979.0  ]

    transpose(x)* inv(A)*x
y[1]^2+y[2]^2+y[3]^2

L2 = zeros(3,3)
n=3
for i in 1:n
    for j in 1:i
        sum1 = 0;
        # summation for diagonals
        if (j == i)
            for k in 1:j
                sum1 += (L2[j,k])^2
                println( """sumAA$(i)$(j) += (L2[$(j),$(k)])^2 """)
            end#for
            L2[j,j] = sqrt(A[j,j] - sum1)
            println(""" L2[$(j),$(j)] = sqrt(A[$(j),$(j)] - sumAA$(i)$(j))   """)

        else
            # Evaluating L(i, j)  using L(j, j)
            for k in 1:j
                sum1 += (L2[i,k] * L2[j,k]);
                println("""sumAA$(i)$(j) += (L2[$(i),$(k)] * L2[$(j),$(k)]) """)
            end#for
            if(L2[j,j] > 0)
                L2[i,j] = (A[i,j] - sum1) /L2[j,j]
                println("""L2[$(i),$(j)] = (A[$(i),$(j)] -sumAA$(i)$(j))/L2[$(j),$(j)] """)

            end#if
        end#if else
    end#for
end#for
L2

#######working unrolled
L2 = zeros(Float32,3,3)
sumAA11=0.0
sumAA21=0.0
sumAA22=0.0
sumAA31=0.0
sumAA32=0.0
sumAA33=0.0

#unrolled 3 by 3 cholesky decomposition# aa
sumAA11 += (L2[1,1])^2
 L2[1,1] = sqrt(A[1,1] - sumAA11)
sumAA21 += (L2[2,1] * L2[1,1])
L2[2,1] = (A[2,1] -sumAA21)/L2[1,1]
sumAA22 += (L2[2,1])^2
sumAA22 += (L2[2,2])^2
 L2[2,2] = sqrt(A[2,2] - sumAA22)
sumAA31 += (L2[3,1] * L2[1,1])
L2[3,1] = (A[3,1] -sumAA31)/L2[1,1]
sumAA32 += (L2[3,1] * L2[2,1])
sumAA32 += (L2[3,2] * L2[2,2])
L2[3,2] = (A[3,2] -sumAA32)/L2[2,2]
sumAA33 += (L2[3,1])^2
sumAA33 += (L2[3,2])^2
sumAA33 += (L2[3,3])^2
 L2[3,3] = sqrt(A[3,3] - sumAA33)

#unrolled forward substitiution
y[1]= x[1]/L2[1,1]
y[2] = (x[2]-L2[2,1]*y[1])/L2[2,2]
y[3]= (x[3]-y[2]*L2[3,2]-y[1]* L2[3,1])/L2[3,3]
#taking square euclidean distance
y[1]^2+y[2]^2+y[3]^2#should be 4.28

#####

L2 = zeros(Float32,3,3)
sumAA11=0.0
sumAA21=0.0
sumAA22=0.0
sumAA31=0.0
sumAA32=0.0
sumAA33=0.0

#unrolled 3 by 3 cholesky decomposition# aa
sumAA11 += (L2[1,1])^2
 L2[1,1] = sqrt(A[1,1] - sumAA11)
sumAA21 += (L2[2,1] * L2[1,1])
L2[2,1] = (A[2,1] -sumAA21)/L2[1,1]
sumAA22 += (L2[2,1])^2
sumAA22 += (L2[2,2])^2
 L2[2,2] = sqrt(A[2,2] - sumAA22)
sumAA31 += (L2[3,1] * L2[1,1])
L2[3,1] = (A[3,1] -sumAA31)/L2[1,1]
sumAA32 += (L2[3,1] * L2[2,1])
sumAA32 += (L2[3,2] * L2[2,2])
L2[3,2] = (A[3,2] -sumAA32)/L2[2,2]
sumAA33 += (L2[3,1])^2
sumAA33 += (L2[3,2])^2
sumAA33 += (L2[3,3])^2
 L2[3,3] = sqrt(A[3,3] - sumAA33)

#unrolled forward substitiution
y[1]= x[1]/L2[1,1]
y[2] = (x[2]-L2[2,1]*y[1])/L2[2,2]
y[3]= (x[3]-y[2]*L2[3,2]-y[1]* L2[3,1])/L2[3,3]
#taking square euclidean distance
y[1]^2+y[2]^2+y[3]^2#should be 4.28

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

L2 = zeros(Float32,3,3)
sumAA11=0.0
sumAA21=0.0
sumAA22=0.0
sumAA31=0.0
sumAA32=0.0
sumAA33=0.0

varianceX= A[1,1]
covarianceXY= A[2,1]
covarianceXZ = A[3,1]
varianceY = A[2,2]
covarianceYZ = A[3,2]
varianceZ = A[3,3]
#unrolled 3 by 3 cholesky decomposition# aa
sumAA11 += (L2[1,1])^2
 L2[1,1] = sqrt(varianceX - sumAA11)
sumAA21 += (L2[2,1] * L2[1,1])
L2[2,1] = (covarianceXY -sumAA21)/L2[1,1]
sumAA22 += (L2[2,1])^2
sumAA22 += (L2[2,2])^2
 L2[2,2] = sqrt(varianceY - sumAA22)
sumAA31 += (L2[3,1] * L2[1,1])
L2[3,1] = (covarianceXZ -sumAA31)/L2[1,1]
sumAA32 += (L2[3,1] * L2[2,1])
sumAA32 += (L2[3,2] * L2[2,2])
L2[3,2] = (covarianceYZ -sumAA32)/L2[2,2]
sumAA33 += (L2[3,1])^2
sumAA33 += (L2[3,2])^2
sumAA33 += (L2[3,3])^2
 L2[3,3] = sqrt(varianceZ - sumAA33)

#unrolled forward substitiution
y[1]= x[1]/L2[1,1]
y[2] = (x[2]-L2[2,1]*y[1])/L2[2,2]
y[3]= (x[3]-y[2]*L2[3,2]-y[1]* L2[3,1])/L2[3,3]
#taking square euclidean distance
y[1]^2+y[2]^2+y[3]^2#should be 4.28

x 1  y 3 z 1   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 0 offset 0
x 2  y 3 z 1   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 0 offset 0
x 3  y 3 z 1   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 0 offset 0
x 1  y 4 z 1   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 0 offset 0
x 2  y 4 z 1   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 0 offset 0
x 3  y 4 z 1   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 0 offset 0
x 1  y 1 z 1   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 0 offset 0
x 2  y 1 z 1   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 0 offset 0
x 3  y 1 z 1   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 0 offset 0
x 1  y 2 z 1   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 0 offset 0
x 2  y 2 z 1   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 0 offset 0
x 3  y 2 z 1   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 0 offset 0
x 1  y 7 z 1   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 0 offset 0
x 2  y 7 z 1   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 0 offset 0
x 3  y 7 z 1   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 0 offset 0
x 1  y 4 z 5   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 0 offset 4
x 2  y 4 z 5   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 0 offset 4
x 3  y 4 z 5   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 0 offset 4
x 1  y 5 z 1   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 0 offset 0
x 2  y 5 z 1   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 0 offset 0
x 3  y 5 z 1   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 0 offset 0
x 1  y 6 z 1   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 0 offset 0
x 2  y 6 z 1   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 0 offset 0
x 3  y 6 z 1   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 0 offset 0
x 1  y 3 z 3   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 0 offset 2
x 2  y 3 z 3   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 0 offset 2
x 3  y 3 z 3   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 0 offset 2
x 1  y 4 z 3   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 0 offset 2
x 2  y 4 z 3   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 0 offset 2
x 3  y 4 z 3   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 0 offset 2
x 1  y 1 z 3   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 0 offset 2
x 2  y 1 z 3   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 0 offset 2
x 3  y 1 z 3   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 0 offset 2
x 1  y 2 z 3   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 0 offset 2
x 2  y 2 z 3   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 0 offset 2
x 3  y 2 z 3   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 0 offset 2
x 1  y 5 z 2   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 1 offset 0
x 2  y 5 z 2   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 1 offset 0
x 3  y 5 z 2   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 1 offset 0
x 1  y 4 z 6   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 1 offset 4
x 2  y 4 z 6   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 1 offset 4
x 3  y 4 z 6   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 1 offset 4
x 1  y 3 z 2   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 1 offset 0
x 2  y 3 z 2   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 1 offset 0
x 3  y 3 z 2   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 1 offset 0
x 1  y 6 z 2   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 1 offset 0
x 2  y 6 z 2   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 1 offset 0
x 3  y 6 z 2   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 1 offset 0
x 1  y 7 z 2   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 1 offset 0
x 2  y 7 z 2   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 1 offset 0
x 3  y 7 z 2   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 1 offset 0
x 1  y 4 z 2   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 1 offset 0
x 2  y 4 z 2   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 1 offset 0
x 3  y 4 z 2   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 1 offset 0
x 1  y 1 z 2   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 1 offset 0
x 2  y 1 z 2   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 1 offset 0
x 3  y 1 z 2   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 1 offset 0
x 1  y 2 z 2   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 1 offset 0
x 2  y 2 z 2   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 1 offset 0
x 3  y 2 z 2   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 1 offset 0
x 1  y 3 z 4   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 1 offset 2
x 2  y 3 z 4   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 1 offset 2
x 3  y 3 z 4   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 1 offset 2
x 1  y 4 z 4   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 1 offset 2
x 2  y 4 z 4   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 1 offset 2
x 3  y 4 z 4   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 1 offset 2
x 1  y 2 z 4   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 1 offset 2
x 2  y 2 z 4   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 1 offset 2
x 3  y 2 z 4   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 1 offset 2
x 1  y 1 z 4   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 1 offset 2
x 2  y 1 z 4   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 1 offset 2
x 3  y 1 z 4   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 1 offset 2
x 4  y 5 z 1   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 0 offset 0
x 5  y 5 z 1   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 0 offset 0
x 6  y 5 z 1   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 0 offset 0
x 4  y 3 z 1   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 0 offset 0
x 5  y 3 z 1   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 0 offset 0
x 6  y 3 z 1   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 0 offset 0
x 4  y 6 z 1   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 0 offset 0
x 5  y 6 z 1   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 0 offset 0
x 6  y 6 z 1   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 0 offset 0
x 4  y 7 z 1   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 0 offset 0
x 5  y 7 z 1   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 0 offset 0
x 6  y 7 z 1   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 0 offset 0
x 4  y 4 z 1   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 0 offset 0
x 5  y 4 z 1   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 0 offset 0
x 6  y 4 z 1   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 0 offset 0
x 4  y 1 z 1   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 0 offset 0
x 5  y 1 z 1   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 0 offset 0
x 6  y 1 z 1   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 0 offset 0
x 4  y 2 z 1   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 0 offset 0
x 5  y 2 z 1   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 0 offset 0
x 6  y 2 z 1   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 0 offset 0
x 4  y 4 z 3   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 0 offset 2
x 5  y 4 z 3   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 0 offset 2
x 6  y 4 z 3   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 0 offset 2
x 4  y 4 z 5   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 0 offset 4
x 5  y 4 z 5   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 0 offset 4
x 6  y 4 z 5   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 0 offset 4
x 4  y 2 z 3   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 0 offset 2
x 5  y 2 z 3   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 0 offset 2
x 6  y 2 z 3   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 0 offset 2
x 4  y 3 z 3   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 0 offset 2
x 5  y 3 z 3   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 0 offset 2
x 6  y 3 z 3   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 0 offset 2
x 4  y 1 z 3   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 0 offset 2
x 5  y 1 z 3   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 0 offset 2
x 6  y 1 z 3   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 0 offset 2
x 4  y 5 z 2   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 1 offset 0
x 5  y 5 z 2   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 1 offset 0
x 6  y 5 z 2   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 1 offset 0
x 4  y 4 z 4   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 1 offset 2
x 5  y 4 z 4   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 1 offset 2
x 6  y 4 z 4   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 1 offset 2
x 4  y 3 z 2   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 1 offset 0
4 5si  y 3 z 2   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 1 offset 0
nx 6  y 3 z 2   xMeta 1  xdim 0 offset 3  idX 3g  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 1 offset 0
x 4  y 6 z 2   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1l  ydim 0 offset 4  idy 2  zMeta 0  zdim 1 offset 0
x 5  y 6e z 2   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 1 offset 0
Vx 6  y 6 z 2   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 1 offset 0
ax 4  y 7 z 2   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 1 offset 0
lx 5  y 7 z 2   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 1 offset 0
[x 6  y 7 z 2   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 1 offset 0
1x 4  y 1 z 4   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 1 offset 2
]x 5  y 1 z 4   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 1 offset 2
)x 6  y 1 z 4   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 1 offset 2
julia> Int64(singleVal[1])x 4  y 1 z 2   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 1 offset 0
x 5  y 1
 z 2   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 1 offset 0
x 6  y 1 z 2   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 1 offset 0
x 4  y 4 z 2   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 1 offset 0
x 5  y 4 z 2   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 1 offset 0
x 6  y 4 z 2   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 1 offset 0
x 4  y 2 z 4   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 1 offset 2
x 5  y 2 z 4   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 1 offset 2
x 6  y 2 z 4   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 1 offset 2
x 4  y 2 z 2   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 1 offset 0
x 5  y 2 z 2   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 1 offset 0
x 6  y 2 z 2   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 1 offset 0
x 4  y 4 z 6   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 1 offset 4
x 5  y 4 z 6   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 1 offset 4
x 6  y 4 z 6   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 1 offset 4
x 4  y 3 z 4   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 1 offset 2
x 5  y 3 z 4   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 1 offset 2
x 6  y 3 z 4   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 1 offset 2
x 7  y 4 z 3   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 0 offset 2
x 7  y 3 z 1   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 0 offset 0
343x 7
  y 5 z 1
   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 0 offset 0
julia> x 7  y 4 z 1   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 0 offset 0
x 7  y 1 z 1   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 0 offset 0
x 7  y 2 z 1   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 0 offset 0
x 7  y 6 z 1   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 0 offset 0
x 7  y 7 z 1   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 0 offset 0
x 7  y 1 z 3   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 0 offset 2
x 7  y 4 z 5   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 0 offset 4
x 7  y 2 z 3   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 0 offset 2
x 7  y 3 z 3   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 0 offset 2
x 7  y 3 z 2   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 0  zdim 1 offset 0
x 7  y 5 z 2   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 0  zdim 1 offset 0
x 7  y 4 z 2   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 0  zdim 1 offset 0
x 7  y 1 z 2   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 0  zdim 1 offset 0
x 7  y 4 z 4   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 1  zdim 1 offset 2
x 7  y 6 z 2   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 0  zdim 1 offset 0
x 7  y 2 z 2   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 0  zdim 1 offset 0
x 7  y 7 z 2   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 0  zdim 1 offset 0
x 7  y 1 z 4   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 1  zdim 1 offset 2
julia> x 7  y 2 z 4   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 1  zdim 1 offset 2
x 7  y 4 z 6   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 2  zdim 1 offset 4
x 7  y 3 z 4   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 1  zdim 1 offset 2
x 1  y 7 z 3   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 0 offset 2
x 2  y 7 z 3   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 0 offset 2
x 3  y 7 z 3   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 0 offset 2
x 1  y 1 z 5   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 0 offset 4
x 2  y 1 z 5   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 0 offset 4
x 3  y 1 z 5   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 0 offset 4
x 1  y 5 z 3   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 0 offset 2
x 2  y 5 z 3   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 0 offset 2
x 3  y 5 z 3   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 0 offset 2
x 1  y 2 z 5   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 0 offset 4
x 2  y 2 z 5   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 0 offset 4
x 3  y 2 z 5   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 0 offset 4
x 1  y 3 z 5   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 0 offset 4
x 2  y 3 z 5   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 0 offset 4
x 3  y 3 z 5   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 0 offset 4
x 1  y 6 z 3   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 0 offset 2
x 2  y 6 z 3   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 0 offset 2
x 3  y 6 z 3   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 0 offset 2
x 1  y 4 z 7   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 3  zdim 0 offset 6
x 2  y 4 z 7   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 3  zdim 0 offset 6
x 3  y 4 z 7   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 3  zdim 0 offset 6
x 1  y 5 z 5   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 0 offset 4
x 2  y 5 z 5   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 0 offset 4
x 3  y 5 z 5   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 0 offset 4
x 1  y 6 z 5   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 0 offset 4
x 2  y 6 z 5   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 0 offset 4
x 3  y 6 z 5   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 0 offset 4
x 1  y 7 z 5   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 0 offset 4
x 2  y 7 z 5   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 0 offset 4
x 3  y 7 z 5   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 0 offset 4
x 1  y 7 z 4   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 1 offset 2
x 2  y 7 z 4   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 1 offset 2
x 3  y 7 z 4   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 1 offset 2
x 1  y 1 z 6   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 1 offset 4
x 2  y 1 z 6   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 1 offset 4
x 3  y 1 z 6   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 1 offset 4
x 1  y 5 z 4   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 1 offset 2
x 2  y 5 z 4   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 1 offset 2
x 3  y 5 z 4   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 1 offset 2
x 1  y 6 z 4   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 1 offset 2
x 2  y 6 z 4   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 1 offset 2
x 3  y 6 z 4   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 1 offset 2
x 4  y 4 z 7   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 3  zdim 0 offset 6
x 5  y 4 z 7   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 4  zMeta 3  zdim 0 offset 6
x 6  y 4 z 7   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 4  zMeta 3  zdim 0 offset 6
x 1  y 2 z 6   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 1 offset 4
x 2  y 2 z 6   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 1 offset 4
x 3  y 2 z 6   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 1 offset 4
x 1  y 5 z 6   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 1 offset 4
x 2  y 5 z 6   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 1 offset 4
x 3  y 5 z 6   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 1 offset 4
x 1  y 3 z 6   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 1 offset 4
x 2  y 3 z 6   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 1 offset 4
x 3  y 3 z 6   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 1 offset 4
x 1  y 6 z 6   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 1 offset 4
x 2  y 6 z 6   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 1 offset 4
x 3  y 6 z 6   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 1 offset 4
x 1  y 7 z 6   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 1 offset 4
x 2  y 7 z 6   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 1 offset 4
x 3  y 7 z 6   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 1 offset 4
x 4  y 7 z 3   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 0 offset 2
x 5  y 7 z 3   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 0 offset 2
x 6  y 7 z 3   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 0 offset 2
x 4  y 1 z 5   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 0 offset 4
x 5  y 1 z 5   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 0 offset 4
x 6  y 1 z 5   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 0 offset 4
x 4  y 5 z 3   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 0 offset 2
x 5  y 5 z 3   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 0 offset 2
x 6  y 5 z 3   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 0 offset 2
x 7  y 4 z 7   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 4  zMeta 3  zdim 0 offset 6
x 4  y 6 z 3   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 0 offset 2
x 5  y 6 z 3   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 0 offset 2
x 6  y 6 z 3   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 0 offset 2
x 4  y 2 z 5   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 0 offset 4
x 5  y 2 z 5   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 0 offset 4
x 6  y 2 z 5   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 0 offset 4
x 4  y 5 z 5   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 0 offset 4
x 5  y 5 z 5   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 0 offset 4
x 6  y 5 z 5   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 0 offset 4
x 4  y 3 z 5   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 0 offset 4
x 5  y 3 z 5   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 0 offset 4
x 6  y 3 z 5   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 0 offset 4
x 4  y 6 z 5   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 0 offset 4
x 5  y 6 z 5   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 0 offset 4
x 6  y 6 z 5   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 0 offset 4
x 4  y 7 z 5   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 0 offset 4
x 5  y 7 z 5   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 0 offset 4
x 6  y 7 z 5   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 0 offset 4
x 4  y 7 z 4   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 1 offset 2
x 5  y 7 z 4   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 1 offset 2
x 6  y 7 z 4   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 1 offset 2
x 4  y 1 z 6   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 1 offset 4
x 5  y 1 z 6   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 1 offset 4
x 6  y 1 z 6   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 1 offset 4
x 4  y 5 z 4   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 1 offset 2
x 5  y 5 z 4   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 1 offset 2
x 6  y 5 z 4   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 1 offset 2
x 4  y 6 z 4   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 1 offset 2
x 5  y 6 z 4   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 1 offset 2
x 6  y 6 z 4   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 1 offset 2
x 4  y 2 z 6   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 1 offset 4
x 5  y 2 z 6   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 1 offset 4
x 6  y 2 z 6   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 1 offset 4
x 4  y 5 z 6   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 1 offset 4
x 5  y 5 z 6   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 1 offset 4
x 6  y 5 z 6   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 1 offset 4
x 4  y 3 z 6   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 1 offset 4
x 5  y 3 z 6   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 1 offset 4
x 6  y 3 z 6   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 1 offset 4
x 4  y 6 z 6   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 1 offset 4
x 5  y 6 z 6   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 1 offset 4
x 6  y 6 z 6   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 1 offset 4
x 4  y 7 z 6   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 1 offset 4
x 5  y 7 z 6   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 1 offset 4
x 6  y 7 z 6   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 1 offset 4
x 7  y 7 z 3   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 0 offset 2
x 7  y 5 z 3   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 0 offset 2
x 7  y 6 z 3   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 0 offset 2
x 7  y 1 z 5   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 0 offset 4
x 7  y 2 z 5   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 0 offset 4
x 7  y 3 z 5   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 0 offset 4
x 7  y 5 z 5   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 0 offset 4
x 7  y 6 z 5   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 0 offset 4
x 7  y 7 z 5   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 0 offset 4
x 7  y 7 z 4   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 1  zdim 1 offset 2
x 7  y 5 z 4   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 1  zdim 1 offset 2
x 7  y 6 z 4   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 1  zdim 1 offset 2
x 7  y 1 z 6   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 2  zdim 1 offset 4
x 7  y 2 z 6   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 2  zdim 1 offset 4
x 7  y 3 z 6   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 2  zdim 1 offset 4
x 7  y 5 z 6   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 2  zdim 1 offset 4
x 7  y 6 z 6   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 2  zdim 1 offset 4
x 7  y 7 z 6   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 2  zdim 1 offset 4
x 1  y 3 z 7   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 3  zdim 0 offset 6
x 2  y 3 z 7   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 3  zdim 0 offset 6
x 3  y 3 z 7   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 3  zdim 0 offset 6
x 1  y 1 z 7   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 3  zdim 0 offset 6
x 2  y 1 z 7   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 3  zdim 0 offset 6
x 3  y 1 z 7   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 3  zdim 0 offset 6
x 1  y 2 z 7   xMeta 0  xdim 0 offset 0  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 3  zdim 0 offset 6
x 2  y 2 z 7   xMeta 0  xdim 0 offset 0  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 3  zdim 0 offset 6
x 3  y 2 z 7   xMeta 0  xdim 0 offset 0  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 3  zdim 0 offset 6
x 1  y 5 z 7   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 3  zdim 0 offset 6
x 2  y 5 z 7   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 3  zdim 0 offset 6
x 3  y 5 z 7   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 3  zdim 0 offset 6
x 1  y 6 z 7   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 3  zdim 0 offset 6
x 2  y 6 z 7   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 3  zdim 0 offset 6
x 3  y 6 z 7   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 3  zdim 0 offset 6
x 1  y 7 z 7   xMeta 0  xdim 0 offset 0  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 3  zdim 0 offset 6
x 2  y 7 z 7   xMeta 0  xdim 0 offset 0  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 3  zdim 0 offset 6
x 3  y 7 z 7   xMeta 0  xdim 0 offset 0  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 3  zdim 0 offset 6
x 4  y 3 z 7   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 3  zdim 0 offset 6
x 5  y 3 z 7   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 3  zMeta 3  zdim 0 offset 6
x 6  y 3 z 7   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 3  zMeta 3  zdim 0 offset 6
x 4  y 1 z 7   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 3  zdim 0 offset 6
x 5  y 1 z 7   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 1  zMeta 3  zdim 0 offset 6
x 6  y 1 z 7   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 1  zMeta 3  zdim 0 offset 6
x 4  y 2 z 7   xMeta 1  xdim 0 offset 3  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 3  zdim 0 offset 6
x 5  y 2 z 7   xMeta 1  xdim 0 offset 3  idX 2  yMeta 0  ydim 0 offset 0  idy 2  zMeta 3  zdim 0 offset 6
x 6  y 2 z 7   xMeta 1  xdim 0 offset 3  idX 3  yMeta 0  ydim 0 offset 0  idy 2  zMeta 3  zdim 0 offset 6
x 4  y 5 z 7   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 3  zdim 0 offset 6
x 5  y 5 z 7   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 1  zMeta 3  zdim 0 offset 6
x 6  y 5 z 7   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 1  zMeta 3  zdim 0 offset 6
x 4  y 6 z 7   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 3  zdim 0 offset 6
x 5  y 6 z 7   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 2  zMeta 3  zdim 0 offset 6
x 6  y 6 z 7   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 2  zMeta 3  zdim 0 offset 6
x 4  y 7 z 7   xMeta 1  xdim 0 offset 3  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 3  zdim 0 offset 6
x 5  y 7 z 7   xMeta 1  xdim 0 offset 3  idX 2  yMeta 1  ydim 0 offset 4  idy 3  zMeta 3  zdim 0 offset 6
x 6  y 7 z 7   xMeta 1  xdim 0 offset 3  idX 3  yMeta 1  ydim 0 offset 4  idy 3  zMeta 3  zdim 0 offset 6
x 7  y 3 z 7   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 3  zMeta 3  zdim 0 offset 6
x 7  y 1 z 7   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 1  zMeta 3  zdim 0 offset 6
x 7  y 2 z 7   xMeta 2  xdim 0 offset 6  idX 1  yMeta 0  ydim 0 offset 0  idy 2  zMeta 3  zdim 0 offset 6
x 7  y 5 z 7   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 1  zMeta 3  zdim 0 offset 6
x 7  y 6 z 7   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 2  zMeta 3  zdim 0 offset 6
x 7  y 7 z 7   xMeta 2  xdim 0 offset 6  idX 1  yMeta 1  ydim 0 offset 4  idy 3  zMeta 3  zdim 0 offset 6


text:
## There are many Poly class methods
## Instances of Poly can be created via:
## @syms x y z
## p = poly(x^2 -x - 1)
##

const SymPoly = Sym

Base.divrem(p::Sym, q::Sym) = _sympy_meth(:div, p, q)
Base.div(p::Sym, q::Sym) = divrem(p,q)[1]
Base.rem(p::Sym, q::Sym) = divrem(p,q)[2]

## rename these, their use is special to polynomials, so we prefix
## Renamed polynomial methods.

## The functions `div`, `rem, `divrem` had slightly different meaning in
## Julia than SymPy, where they are polynomial methods. Here we rename
## them by prefixing with "poly". We do the same for `roots`, so as to
## not conflict with the `roots` function from `Polynomials.jl`.

## [This isn't really needed. We do keep the `polyroots`, as `roots` does not live in base.]

"""

Polynomial division. Renamed from `div` in SymPy to avoid confusion with Julia's `div`. Also, same as `divrem`.

"""
polydiv(ex::Sym, args...; kwargs...) = sympy_meth(:div, ex, args...; kwargs...)

"""

Polynomial division remainer.

"""
polyrem(ex::Sym, args...; kwargs...) = rem(ex, args...; kwargs...)

"""

Polynomial division with remainder. Renamed from `divrem` in SymPy to avoid confusion with Julia's `divrem`

"""
polydivrem(ex::Sym, args...; kwargs...) = _sympy_meth(:div, ex, args...; kwargs...)

"""

Find roots of a polynomial. Renamed from `roots` in
SymPy to avoid confusion with the `roots` function of `Polynomials`

"""
polyroots(ex::Sym, args...; kwargs...) = sympy_meth(:roots, ex, args...; kwargs...)
export polydiv, polyrem, polydivrem, polyroots

polynomial_sympy_methods = (
                            :groebner,
                            :solve_poly_system,
                            :resultant,
                            :cancel,
                            :apart, :together,
                            :poly,
                            :poly_from_expr,
                            :degree_list,
                            :pdiv, :prem, :pquo, :pexquo,
                            :quo, :exquo,
                            :half_gcdex, :gcdex,
                            :invert,
                            :subresultants,
                            :discriminant,
                            :terms_gcd,
                            :cofactors,
                            :gcd_list,
                            :lcm_list,
                            :monic,
                            :content,
                            :compose, :decompose,
                            :sturm,
                            :gff_list, :gff,
                            :sqf_norm, :sqf_part, :sqf_list, :sqf,
                            :factor_list,
                            :intervals,
                            :refine_root,
                            :count_roots,
                            :real_roots,
                            :nroots,
                            :ground_roots,
                            :reduced,
                            :is_zero_dimensional,
                            :symmetrize,
                            :horner,
                            :viete,
                            :construct_domain,
                            :minimal_polynomial,
                            :minpoly, :primitive_element,
                            :field_isomorphism, :to_number_field
#                            :numer#, :denom
#                            :roots ## conflict with Roots.roots and functionality provided by solve
                            )

degree(x::Sym;kwargs...) = _sympy_meth(:degree, x; kwargs...)
export degree

denom(x::Sym) = _sympy_meth(:denom, x)
numer(x::Sym) = _sympy_meth(:numer, x)
export denom, numer
Base.denominator(x::SymbolicObject) = denom(x)
Base.numerator(x::SymbolicObject) = numer(x)

polynomial_sympy_methods_base = (:expand,
                                 :factor #,
                                    #:trunc
                                    )

## SymPy Poly class

"""

Constructor for polynomials. This allows certain polynomial-specifec
methods to be called, such as `coeffs`. There is also `poly` that does
the same thing. (We keep `Poly`, as it is named like a constructor.)

Examples:
```
p = Poly(x^2 + 2x + 1, x)  #  a poly in x over Z
coeffs(p)  # [1, 2, 1]
p = Poly(a*x^2 + b*x + c, x)  # must specify variable that polynomial is over.
coeffs(p)  ## [a,b,c]
```

"""
Poly(x::Sym) = sympy."Poly"(x)
Poly(x::Sym, args...; kwargs...) = sympy_meth(:Poly, x, args...; kwargs...)
export Poly

## Poly class methods that aren't sympy methods
## e.g., must be called as obj.method(...), but not method(obj, ...)
## not coeffs(x^2 -1), rather poly(x^2 - 1)  |> coeffs
## or poly(a*x^2 + b*x + c, x) |> coeffs to specify a specific variable
polynomial_instance_methods = (:EC, :ET, :LC, :LM, :LT, :TC, ## no .abs()
                               :add_ground,
                               :all_coeffs,
                              :all_monoms, :all_roots, :all_terms,
                              :as_dict, :as_expr, :as_list,
                              :clear_denoms,
                              :deflate,
                              :exclude,
                              :l1_norm,
                              :lift,
                              :ltrim,
                              :max_norm,
                              :monoms,
                              :mul_ground,
                              :neg,
                              :nth, :nth_power_roots_poly,
                              :per,
                              :rat_clear_denoms,
                              :retract,
                              :root,
                              :set_domain, :set_modulus,
                              :sqr, :sub_ground,
                              :to_exact, :to_field, :to_ring, :total_degree,
                              :unify
)

"""

    coeffs(p, [x])

Return coefficients, `[a0, a1, ..., an]`, of a polynomial `p = a0 + a1*x + ... + an*x^n` in `D[x]`.

Note: SymPy has the `all_coeffs` function to do this for objects of class `Poly`. The SymPy `coeffs` function returns the non-zero values. Use `sympy"coeffs"` for this functionality.

This function returns the coefficients starting with the constant and
includes zeros, when present. This matches the order of the `coeffs`
function of the `Polynomials` package. This is implemented as
`reverse(all_coeffs(Poly(p,[x])))`.

The expression is converted to a object of class `Poly`. This requires a choices of symbol after which the expression is viewed as a univariate polynomial over a domain. If no symbol is passed in, a heuristic is used, for convenience: when only one free symbol is present it will be used; if a symbol "x" is present it will be used; otherwise the first symbol returned by `free_symbols` will be used. This will clearly fail for constant polynomials.

Examples:
```
@vars x a b c
ex = 2x^2 + 3x + 4
coeffs(ex)  # [4, 3, 2], same as `coeffs(ex, x)`

ex = a*x^2 + b*x + c
coeffs(ex)  # [c, b, a], same as `coeffs(ex, x)`

## may need to be explicit if `x` not present:
@vars y
ex = a*y^2 + b*y + c
coeffs(ex)  # not [c,b,a] as `c` is first free symbol
coeffs(ex, y)  # [c, b, a]
```
"""
function coeffs(p::Sym, x::Sym)
    q = Poly(p, x)
    reverse(object_meth(q, :all_coeffs))
end
function coeffs(p::Sym)
    vars = free_symbols(p)
    length(vars) == 0 && DomainError("No free variables specified")
    if length(vars) == 1
        return coeffs(p, vars[1])
    else
        for var in vars
            if "x" == string(var)
                return coeffs(p, var)
            end
        end
        return coeffs(p, vars[1])
    end
end
export coeffs

polynomial_predicates = (
                         :is_cyclotomic,
                         :is_ground,
                         :is_homogeneous,
                         :is_irreducible,
                         :is_linear,
                         :is_monic,
                         :is_monomial,
                         :is_multivariate,
                         :is_one,
#                        :is_primitive,
                         :is_quadratic,
                         :is_sqf,
                         :is_univariate,
                         :is_zero)

is_primitive(x::Sym) = x.is_primitive
export is_primitive

##
"""

    interpolate([xs], ys, x)

Find interpolating polynomial in the variable x through the points (xs[1], ys[1]), ..., (xs[n], ys[n])

Examples
```jltest
@vars x
interpolate(sin.(1:4), x)          # uses xs = 1:length(ys(
interpolate(1:4, sin.(1:4), x)     # xs, ys
interpolate([(1,2), (2,3), (3,2)], x)  # [(x1,y1), ..., (xn, yn)]
```
"""
interpolate(xsys::AbstractVector, x) = _sympy_meth(:interpolate, xsys, x)
interpolate(xs::AbstractVector, ys::AbstractVector, x::Sym) = interpolate(collect(zip(xs, ys)), x)
export interpolate


text:
# This file was generated, do not modify it. # hide
@pipeline LinMod(fs = FeatureSelector(features=[:hp1]),
                 lr = LinearRegressor());


text:
function fRhoBGrd(x,Param)
if lower(Param.ProfRhoBGrd)== "baldaufcart"
    delta=Param.Grav/(Param.Rd*Param.T0);
    p=Param.p0*exp(-delta*x[3]);
    TLoc=Param.T0;
    Rho=p/(Param.Rd*TLoc);
else
    Rho=0;
end
return Rho
end


text:
# This file was generated by the Julia Swagger Code Generator
# Do not modify this file directly. Modify the swagger specification instead.

struct ExpressRouteCircuitPeeringsApi <: SwaggerApi
    client::Swagger.Client
end

"""

Creates or updates a peering in the specified express route circuits.
Param: resourceGroupName::String (required)
Param: circuitName::String (required)
Param: peeringName::String (required)
Param: peeringParameters::ExpressRouteCircuitPeering (required)
Param: api_version::String (required)
Param: subscriptionId::String (required)
Return: ExpressRouteCircuitPeering
"""
function _swaggerinternal_expressRouteCircuitPeeringsCreateOrUpdate(_api::ExpressRouteCircuitPeeringsApi, resourceGroupName::String, circuitName::String, peeringName::String, peeringParameters, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "PUT", ExpressRouteCircuitPeering, "/subscriptions/{subscriptionId}/resourceGroups/{resourceGroupName}/providers/Microsoft.Network/expressRouteCircuits/{circuitName}/peerings/{peeringName}", ["azure_auth"], peeringParameters)
    Swagger.set_param(_ctx.path, "resourceGroupName", resourceGroupName)  # type String
    Swagger.set_param(_ctx.path, "circuitName", circuitName)  # type String
    Swagger.set_param(_ctx.path, "peeringName", peeringName)  # type String
    Swagger.set_param(_ctx.path, "subscriptionId", subscriptionId)  # type String
    Swagger.set_param(_ctx.query, "api-version", api_version)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/json"] : [_mediaType])
    return _ctx
end

function expressRouteCircuitPeeringsCreateOrUpdate(_api::ExpressRouteCircuitPeeringsApi, resourceGroupName::String, circuitName::String, peeringName::String, peeringParameters, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = _swaggerinternal_expressRouteCircuitPeeringsCreateOrUpdate(_api, resourceGroupName, circuitName, peeringName, peeringParameters, api_version, subscriptionId; _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function expressRouteCircuitPeeringsCreateOrUpdate(_api::ExpressRouteCircuitPeeringsApi, response_stream::Channel, resourceGroupName::String, circuitName::String, peeringName::String, peeringParameters, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = _swaggerinternal_expressRouteCircuitPeeringsCreateOrUpdate(_api, resourceGroupName, circuitName, peeringName, peeringParameters, api_version, subscriptionId; _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

"""

Deletes the specified peering from the specified express route circuit.
Param: resourceGroupName::String (required)
Param: circuitName::String (required)
Param: peeringName::String (required)
Param: api_version::String (required)
Param: subscriptionId::String (required)
Return: Nothing
"""
function _swaggerinternal_expressRouteCircuitPeeringsDelete(_api::ExpressRouteCircuitPeeringsApi, resourceGroupName::String, circuitName::String, peeringName::String, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "DELETE", Nothing, "/subscriptions/{subscriptionId}/resourceGroups/{resourceGroupName}/providers/Microsoft.Network/expressRouteCircuits/{circuitName}/peerings/{peeringName}", ["azure_auth"])
    Swagger.set_param(_ctx.path, "resourceGroupName", resourceGroupName)  # type String
    Swagger.set_param(_ctx.path, "circuitName", circuitName)  # type String
    Swagger.set_param(_ctx.path, "peeringName", peeringName)  # type String
    Swagger.set_param(_ctx.path, "subscriptionId", subscriptionId)  # type String
    Swagger.set_param(_ctx.query, "api-version", api_version)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/json"] : [_mediaType])
    return _ctx
end

function expressRouteCircuitPeeringsDelete(_api::ExpressRouteCircuitPeeringsApi, resourceGroupName::String, circuitName::String, peeringName::String, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = _swaggerinternal_expressRouteCircuitPeeringsDelete(_api, resourceGroupName, circuitName, peeringName, api_version, subscriptionId; _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function expressRouteCircuitPeeringsDelete(_api::ExpressRouteCircuitPeeringsApi, response_stream::Channel, resourceGroupName::String, circuitName::String, peeringName::String, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = _swaggerinternal_expressRouteCircuitPeeringsDelete(_api, resourceGroupName, circuitName, peeringName, api_version, subscriptionId; _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

"""

Gets the specified peering for the express route circuit.
Param: resourceGroupName::String (required)
Param: circuitName::String (required)
Param: peeringName::String (required)
Param: api_version::String (required)
Param: subscriptionId::String (required)
Return: ExpressRouteCircuitPeering
"""
function _swaggerinternal_expressRouteCircuitPeeringsGet(_api::ExpressRouteCircuitPeeringsApi, resourceGroupName::String, circuitName::String, peeringName::String, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "GET", ExpressRouteCircuitPeering, "/subscriptions/{subscriptionId}/resourceGroups/{resourceGroupName}/providers/Microsoft.Network/expressRouteCircuits/{circuitName}/peerings/{peeringName}", ["azure_auth"])
    Swagger.set_param(_ctx.path, "resourceGroupName", resourceGroupName)  # type String
    Swagger.set_param(_ctx.path, "circuitName", circuitName)  # type String
    Swagger.set_param(_ctx.path, "peeringName", peeringName)  # type String
    Swagger.set_param(_ctx.path, "subscriptionId", subscriptionId)  # type String
    Swagger.set_param(_ctx.query, "api-version", api_version)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/json"] : [_mediaType])
    return _ctx
end

function expressRouteCircuitPeeringsGet(_api::ExpressRouteCircuitPeeringsApi, resourceGroupName::String, circuitName::String, peeringName::String, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = _swaggerinternal_expressRouteCircuitPeeringsGet(_api, resourceGroupName, circuitName, peeringName, api_version, subscriptionId; _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function expressRouteCircuitPeeringsGet(_api::ExpressRouteCircuitPeeringsApi, response_stream::Channel, resourceGroupName::String, circuitName::String, peeringName::String, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = _swaggerinternal_expressRouteCircuitPeeringsGet(_api, resourceGroupName, circuitName, peeringName, api_version, subscriptionId; _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

"""

Gets all peerings in a specified express route circuit.
Param: resourceGroupName::String (required)
Param: circuitName::String (required)
Param: api_version::String (required)
Param: subscriptionId::String (required)
Return: ExpressRouteCircuitPeeringListResult
"""
function _swaggerinternal_expressRouteCircuitPeeringsList(_api::ExpressRouteCircuitPeeringsApi, resourceGroupName::String, circuitName::String, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = Swagger.Ctx(_api.client, "GET", ExpressRouteCircuitPeeringListResult, "/subscriptions/{subscriptionId}/resourceGroups/{resourceGroupName}/providers/Microsoft.Network/expressRouteCircuits/{circuitName}/peerings", ["azure_auth"])
    Swagger.set_param(_ctx.path, "resourceGroupName", resourceGroupName)  # type String
    Swagger.set_param(_ctx.path, "circuitName", circuitName)  # type String
    Swagger.set_param(_ctx.path, "subscriptionId", subscriptionId)  # type String
    Swagger.set_param(_ctx.query, "api-version", api_version)  # type String
    Swagger.set_header_accept(_ctx, ["application/json"])
    Swagger.set_header_content_type(_ctx, (_mediaType === nothing) ? ["application/json"] : [_mediaType])
    return _ctx
end

function expressRouteCircuitPeeringsList(_api::ExpressRouteCircuitPeeringsApi, resourceGroupName::String, circuitName::String, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = _swaggerinternal_expressRouteCircuitPeeringsList(_api, resourceGroupName, circuitName, api_version, subscriptionId; _mediaType=_mediaType)
    Swagger.exec(_ctx)
end

function expressRouteCircuitPeeringsList(_api::ExpressRouteCircuitPeeringsApi, response_stream::Channel, resourceGroupName::String, circuitName::String, api_version::String, subscriptionId::String; _mediaType=nothing)
    _ctx = _swaggerinternal_expressRouteCircuitPeeringsList(_api, resourceGroupName, circuitName, api_version, subscriptionId; _mediaType=_mediaType)
    Swagger.exec(_ctx, response_stream)
end

export expressRouteCircuitPeeringsCreateOrUpdate, expressRouteCircuitPeeringsDelete, expressRouteCircuitPeeringsGet, expressRouteCircuitPeeringsList


text:
### A Pluto.jl notebook ###
# v0.12.10

using Markdown
using InteractiveUtils

# ╔═╡ e077ab16-22e5-11eb-15ed-5f0654bf7328
using Pkg, DrWatson

# ╔═╡ e75a9a06-22e5-11eb-2e36-8d4f62b830ed
begin
	@quickactivate "StatisticsWithJuliaPlutoNotebooks"
	using Random, Distributions, StatsPlots
	Random.seed!(0)
end

# ╔═╡ ca5fbdc8-22e5-11eb-0a60-cb0a127a3ca5
md"## Listing4.15"

# ╔═╡ fa68607e-22e5-11eb-0558-c9a4d9f77426
begin
	mu1, sigma1 = 10, 5
	mu2, sigma2 = 40, 12
	dist1, dist2 = Normal(mu1,sigma1), Normal(mu2,sigma2)
	p = 0.3
	mixRv() = (rand() <= p) ? rand(dist1) : rand(dist2)

	n = 2000
	data = [mixRv() for _ in 1:n]
end

# ╔═╡ cd520f9c-2509-11eb-36ce-8340182a95eb
begin
	density(data, c=:blue, label="Density via StatsPlots",
			xlims=(-20,80), ylims=(0,0.035))
	stephist!(data, bins=50, c=:black, norm=true,
		label="Histogram", xlabel="x", ylabel = "Density")
end

# ╔═╡ 475ff888-22e6-11eb-2354-09f8f40a8e12
md"## End of listing4.15"

# ╔═╡ Cell order:
# ╟─ca5fbdc8-22e5-11eb-0a60-cb0a127a3ca5
# ╠═e077ab16-22e5-11eb-15ed-5f0654bf7328
# ╠═e75a9a06-22e5-11eb-2e36-8d4f62b830ed
# ╠═fa68607e-22e5-11eb-0558-c9a4d9f77426
# ╠═cd520f9c-2509-11eb-36ce-8340182a95eb
# ╟─475ff888-22e6-11eb-2354-09f8f40a8e12


text:
using Maracas
import Maracas: AbstractTestSet, rm_spec_char, rpad_title
@testset "Maracas" begin
    @testset "'describe' returns a TestSet" begin
        @test isa(@describe("description", begin end), AbstractTestSet)
    end

    @testset "'describe' TestSet contains description" begin
        expected = "expected description"
        ts = @describe("expected description", begin end)
        @test occursin( expected, ts.description)
    end

    @testset "'describe' TestSet description is colored within 'title' env var" begin
        ts = @describe("description", begin end)
        @test occursin( MARACAS_SETTING[:title], ts.description)
    end

    @testset "'describe' TestSet description color can be changed" begin
        title_color = MARACAS_SETTING[:title]
        set_title_style(:blue)
        ts = @describe("description", begin end)
        @test occursin( Base.text_colors[:blue], ts.description)
        MARACAS_SETTING[:title] = title_color
    end

    @testset "'it' returns a TestSet" begin
        @test isa(@it("description", begin end), AbstractTestSet)
    end

    @testset "'it' TestSet contains description" begin
        expected = "expected description for it"
        ts = @it("expected description for it", begin end)
        @test occursin( expected, ts.description)
    end

    @testset "'it' TestSet contains [Spec]" begin
        ts = @it("description", begin end)
        @test occursin( "[Spec]", ts.description)
    end

    @testset "'it' TestSet description is colored within 'spec' env var" begin
        ts = @it("description", begin end)
        @test occursin( MARACAS_SETTING[:spec], ts.description)
    end

    @testset "'test' returns a TestSet" begin
        @test isa(@unit("description", begin end), AbstractTestSet)
    end

    @testset "'test' TestSet contains description" begin
        expected = "expected description for unit"
        ts = @unit("expected description for unit", begin end)
        @test occursin( expected, ts.description)
    end

    @testset "'test' TestSet contains [Test]" begin
        ts = @unit("description", begin end)
        @test occursin( "[Test]", ts.description)
    end

    @testset "'test' TestSet description is colored within 'test' env var" begin
        ts = @unit("description", begin end)
        @test occursin( MARACAS_SETTING[:test], ts.description)
    end

    @testset "padding results and descriptions" begin
        @testset "remove special chars returns empty string when special char is given" begin
            @test rm_spec_char(MARACAS_SETTING[:spec]) == ""
        end
        @testset "remove special chars returns a string with special chars removed" begin
            given = string("maracas : ", MARACAS_SETTING[:spec], "do tchik tchik tchik")
            expected = "maracas : do tchik tchik tchik"
            @test rm_spec_char(given) == expected
        end
        @testset "rpad_title of a special char returns :title_length spaces" begin
            given = MARACAS_SETTING[:test]
            expected = " "^MARACAS_SETTING[:title_length]
            @test rm_spec_char(rpad_title(given)) == expected
        end
        @testset "test title is cut when too long" begin
            given = "-"^(MARACAS_SETTING[:title_length] + 10)
            @test length(rm_spec_char(rpad_title(given))) == MARACAS_SETTING[:title_length]
        end
        @testset "test title end is replaced with ellipsis when too long" begin
            given = "-"^(MARACAS_SETTING[:title_length] + 10)
            @test rpad_title(given)[end-3:end-1] == "..."
        end

    end
end

include("doc_examples.jl")


text:
## ExoplanetsSysSim/src/target.jl
## (c) 2015 Eric B. Ford

#using Distributions

immutable KeplerTarget
  #sys::PlanetarySystem   # Make array for planetary systems aroud multiple stars in one target?
  sys::Vector{PlanetarySystemAbstract}
  cdpp::Array{Float64,2} # fractional, not ppm; 2D to allow for multiple time scales, months, quarters or seasons/spacecraft rotation angles
                           # QUERY: Do we want a separate CDPP for SC?  Or will CDPP's in different months be for LC/SC depending on this variable?
                            # QUERY: Should this be moved to KeplerTargetObs?
  contam::Float64     # QUERY: Do we want/need this, since we're able to generate multiple stars in a single target?
  data_span::Float64
  duty_cycle::Float64
  #channel::Int64         # E.g., if we cared which Kepler channel the target fell on
  #has_sc::Vector{Bool}   # TODO OPT: Make Immutable Vector or BitArray for speed?  QUERY: Should this go in KeplerTargetObs?
  #                       # QUERY: Do we want a separate CDPP for SC?  Or will CDPP's in different months be for LC/SC depending on this variable?
  #ra::Float64           # E.g., if we cared about position on sky     QUERY:  Should we replace with galactic longitude and latitute?
  #dec::Floa64           #
end
num_planets(t::KeplerTarget) = sum( num_planets, t.sys)
flux(t::KeplerTarget) = sum(flux,t.sys)+t.contam

function draw_asymmetric_normal(mu::Real, sig_plus::Real, sig_minus::Real)
  stdn = rand()
  mu + ( (stdn>=0.0) ? sig_plus*stdn : sig_minus*stdn )
end

function generate_kepler_target_from_table(sim_param::SimParam)
  #const  generate_star = get_function(sim_param,"generate_star")
  const  generate_planetary_system = get_function(sim_param,"generate_planetary_system")
  const use_star_table_sigmas = false

  max_star_id = num_usable_in_star_table()
  @assert(1<=max_star_id)
  star_id = rand(1:max_star_id)
  mass = 0.0
  dens = 0.0
  radius = 0.0
  if use_star_table_sigmas
    while !(0.0<mass<100.0)
       mass = draw_asymmetric_normal( star_table(star_id,:mass), star_table(star_id,:mass_err1), star_table(star_id,:mass_err2) )
    end
    while !(0.0<dens<1000.0)
      dens = draw_asymmetric_normal( star_table(star_id,:dens), star_table(star_id,:dens_err1), star_table(star_id,:dens_err2) )
    end
    while !(0.0<radius<10000.0)
      radius = draw_asymmetric_normal( star_table(star_id,:radius), star_table(star_id,:radius_err1), star_table(star_id,:radius_err2) )
    end
  else
    mass   = star_table(star_id,:mass)
    dens   = star_table(star_id,:dens)
    radius = star_table(star_id,:radius)
  end
  star = SingleStar(radius,mass,1.0, star_id)        # TODO SCI: Allow for blends, binaries, etc.
  cdpp = 1.0e-6 * star_table(star_id, :rrmscdpp04p5) * sqrt(4.5/24.0 / LC_duration )  # TODO SCI: Allow for multiple timescales
  contam = 0.0 # rand(LogNormal(1.0e-3,1.0))      # TODO SCI: Come up with better description of Kepler targets, maybe draw from real contaminations
  data_span = star_table(star_id, :dataspan)
  duty_cycle = star_table(star_id, :dutycycle)
  # ch = rand(DiscreteUniform(1,84))
  ps = generate_planetary_system(star, sim_param)
  return KeplerTarget([ps],fill(cdpp,num_cdpp_timescales,num_quarters),contam,data_span,duty_cycle) #,ch )
end

function generate_kepler_target_simple(sim_param::SimParam)
  const  generate_star = get_function(sim_param,"generate_star")
  const  generate_planetary_system = get_function(sim_param,"generate_planetary_system")
  const star::StarAbstract = generate_star(sim_param)

  const mean_log_cdpp = 4.9759601617565465   # mean frmo star table
  const stddev_log_cdpp = 0.6704860437536709 # std dev from star table
  rrmscdpp_5hr = exp(mean_log_cdpp+stddev_log_cdpp*randn())
  cdpp = 1.0e-6 * rrmscdpp_5hr * sqrt(5.0/24.0 / LC_duration )
  contam = 0.0 # rand(LogNormal(1.0e-3,1.0))   # TODO SCI: Come up with better description of Kepler targets, maybe draw from real contaminations
  #ch = rand(DiscreteUniform(1,84))
  ps = generate_planetary_system(star, sim_param)
  return KeplerTarget([ps],fill(cdpp,num_cdpp_timescales,num_quarters),contam,mission_data_span,mission_duty_cycle) #,ch)
end

function test_target(sim_param::SimParam)
  generate_kepler_target_simple(sim_param)
  setup_star_table(sim_param)
  generate_kepler_target_from_table(sim_param)
end


text:
function _validate_children(children::Union{Vector, Tuple}, ids::Set{Symbol})
    for child in children
        _validate(child, ids)
    end
end
function _validate_children(children::Component, ids::Set{Symbol})
    _validate(children, ids)
end
function _validate_children(children, ids::Set{Symbol}) end

function _validate(comp::Component, ids::Set{Symbol})
    if hasproperty(comp, :id) && !isnothing(comp.id)
        id = Symbol(comp.id)
        id in ids && error("Duplicate component id found in the initial layout: $(id)")
        push!(ids, id)
    end
    hasproperty(comp, :children) && _validate_children(comp.children, ids)
end
function _validate(non_comp, ids::Set{Symbol}) end

function validate(comp::Component)
    _validate(comp, Set{Symbol}())
end